STARTUP AND STABILITY OF A GASEOUS CORE
NUCLEAR REACTOR SYSTEM
By
KIRATADAS KUTIKKAD
A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL
FULFILLMENT OF THE REQUIREMENTS FOR
THE DEGREE OF
DOCTOR OF PHILOSOPHY
UNIVERSITY OF FLORIDA
1991
ACKNOWLEDGMENT'S
The author wishes to express his sincere gratitude and
appreciation to Dr. Edward Dugan for his patience. constant guidance and encouragement during the course of this investigation. The author has profited immensely by associating professionally with Dr. Dugan. Dr. Dugan's deep understanding of the physics of nuclear reactors and, in particular, his knowledge of reactor kinetics has been a constant source of help to the author during the course of this work.
The author also wishes to extend his sincere thanks to the
members of his supervisory committee for their willingness to assist in the preparation of this dissertation. Dr. Alan Jacobs needs special mention in this regard. His comments and criticisms were of great help during the final phases of this work.
Financial assistance was provided partly by the Nuclear
Engineering Department and partly by the Innovative Nuclear Space Power Institute at the University of Florida. The author is thankful to Dr. Nils Diaz, Dr. Ed Dugan, and Prof. Jim Tulenko for this financial support.
The author wishes to express his appreciation to a number of his fellow graduate students at the University of Florida. Special thanks go
to Dr. Mathew Panicker and Jerry Welch for the many technical discussions the author had with these two individuals, and to Bob Williams for his cheerful company.
Appreciation and affection are extended to the author's mother, sisters, sistersinlaw, brothers and brothersinlaw for their constant support and encouragement. Deepest appreciation also goes to his friends whose company made the long time it took to finish this work much more enjoyable, quite memorable and precious.
Finally, special thanks are extended to Wale Oladiran, Paul Miller and Charlie McKibben of the University of Missouri Research Reactor for their help and understanding during the final preparation of this manuscript.
TABLE OF CONTENTS
ACKNOWLEDGEMENTS iii
ABSTRACT vii
CHAPTERS
1 INTRODUCTION1I
Advantages of Gaseous Fuel Over Solid Fuel  2
Special Advantages of Gas Core Reactors
for Space Power Applications 3
Present Work Description 5
Dissertation Organization 9
2 PREVIOUS WORK ON GASEOUS CORE REACTORS 11
Previous Gas Core Reactor Work at the
University of Florida 14
Present System Description 16
Summary of Previous Work Done in
Related Areas 20
3 ANALYSIS OF TRANSITION TO HIGH POWER MODE 23
System Description 25
Selection of External Moderator Thickness  27 Spherical Versus Cylindrical Geometry31
Estimation of Gas Core Reactor Reactivity
Coefficients Using the Equivalent Spherical Geometry 36
Models Employed in the Computer Program  45 Programming Aspects 68
Testing of the Computer Program 70
Results 79
Page
4 BPGCR STABILITY ANALYSIS 104
BPGCR Linear Stability Analysis 135
Inherent Stability of the BPGCR Attached to
an ExternalLoop 143
5 CONCLUSIONS, RECOMMENDATIONS AND AREAS
OF FURTHER RESEARCH 171
Results 172
Recommnendations for Future Work 176 APPENDICES
A THE COMPUTER PROGRAM DYNAM 182
B BPGCR SYSTEM EQUATIONS AND THE DERIVATION
OF THE BPGCR TRANSFER FUNCTIONS 209 LIST OF REFERENCES 230
BIOGRAPHICAL SKETCH 233
Abstract of Dissertation Presented to the Graduate
School of the University of Florida in Partial Fulfillment
of the Requirements for the Degree of Doctor of Philosophy
STARTUP AND STABILITY OF A GASEOUS CORE NUCLEAR REACTOR SYSTEM By
Kiratadas Kutikkad
December, 1991
Chairman: Dr. Edward T. Dugan Major Department: Nuclear Engineering Sciences
Dynamic aspects of a conceptual Gas Core Reactor (GCR) are investigated in this dissertation. The reactor is a bimodal, closed cycle, circulating fuel system capable of operating in low and high power modes. It uses uranium hexafluoride (UF6) in gaseous form as the fuel, and the fuel gas mixture of UF6 and helium also serves as the coolant/working fluid.
The gaseous nature of the fuel, plus the fact that this is a
circulating fuel reactor, makes the system unique, and the dynamic behavior tends to be quite different from that of conventional, solid core systems. The system behavior while undergoing rapid transition from a low power mode of operation to a high power mode of
operation is investigated with the help of a computer program developed for that purpose. The program, DYNAM, written in standard FORTRAN 77, simultaneously solves the relevant point reactor kinetics, thermodynamic, heat transfer, and 1D isentropic flow equations of the GCR.
A nonlinear stability analysis of the system is also performed
using DYNAM. The results of this analysis indicate that the investigated GCR is an inherently stable systemmainly due to the strong negative fuel mass feedback effect. Due to this strong fuel mass feedback, the response of the GCR to external reactivity perturbations is found to be quite different from that of conventional, solid core reactors. For example, at full power operation, it is observed that the inherent reactor response is a power decrease for a positive, as well as a negative, external step reactivity insertion. To increase GCR power, it is necessary to provide additional positive reactivity during the transient by varying a system operating parameter (for example, by increasing the core inlet pressure).
A linear stability analysis of the GCR shows that the linear
analysis is often adequate for correctly predicting system stability. However, for certain operating regimes, where the fuel mass feedback coefficient is large, the linear analysis breaks down due to the strong nonlinear effects. For these regimes, the nonlinear analysis has to be performed to reliably predict the system stability.
viii
CHAPTER 1
INTRODUCTION
Terrestrial power generation making use of nuclear fission is currently accomplished with solid fueled reactors. The fuel usually is in the form of small pellets made of uranium dioxide (U02) or in some cases uranium metal or some other compounds of uranium like uranium carbide or nitride. The concept of using nuclear fuel in gaseous form was first suggested in 1955 (Ref. 1). Here the fuel can be either uranium plasma (provided the temperature in the system is kept high enough for the plasma to exist) or a compound like uranium hexafluoride (UF6) or uranium tetrafluoride (UF4) in gaseous form. Nuclear fuel in gaseous form has many potential advantages over that in solid form which could eventually make it the fuel form of choice for many special applications including space power.
In this dissertation, dynamic aspects of a conceptual gas core
reactor system are investigated. The reactor uses highly enriched U176 or UF4 in gaseous form as the fuel. The peak fuel temperatures attained are in the range of only a few thousands of degree Kelvin and, hence, the fuel is not in a plasma state.
Knowledge of a reactor system's behavior in the nonsteady state is of primary importance in the design and construction of any new type of reactor. The purpose of the research described in this dissertation is twofold. The first goal is to perform the system startup analysis, i.e., analysis and prediction of the system behavior during various startup scenarios and the second goal is to analyze the stability of the system's high power mode of operation, i.e., to establish the conditions under which the reactor can be operated in a safe and reliable manner at high power levels.
Advantages of Gaseous Fuel Over Solid Fuel
One of the major limitations of conventional solid core
pressurized water reactors (PWR) or boiling water reactors (BWR) is their relatively low working fluid temperature. Since the zircaloy cladding of the fuel loses its strength at around 700 K, the working fluid/coolant temperature has to be kept below this value. This is a low value compared to the 1500 K or so reached at the centerline of an average fuel pellet inside the core, and this results in a low power plant energy conversion efficiency compared to a coal or oil fired power plant. The reason for not replacing zircaloy with some other metal or alloy with better mechanical strength is that, except for this drawback, zircaloy is an excellent cladding material with many desirable properties like low neutron absorption cross section, good resistance to corrosion by water at elevated temperatures, and reasonably good resistance to radiation damage. By going to a gaseous
3
fuel, this serious limitation of solid core systems can be circumvented. A gaseous core system can operate at extremely high temperatures (and high efficiencies) since the fuel, either by itself or when mixed with another gas, can also be the coolant/working fluid, and, unlike in a solid core system, there may not need to be any heat transfer process from the fuel to a coolant/working fluid. Practical heat removal and materials considerations generally force one to operate the gas core reactor at temperatures in the range of a few thousands of degree Kelvin. An exception to this is the plasma core reactors where the fuel temperatures are in the tens of thousands of degree Kelvin range. Here special arrangements are made to keep the fuel away from the containment wall.
A gaseous core system also opens up the possibility of continuous operation of the power plant with online refuelling and online fission product removal. Other advantages associated with gas core reactors include a much less complicated core structure, high fuel utilization, simple fuel management, and inherent safety due to the negative reactivity associated with an expanding gaseous fuel.
Special Advantages of Gas Core Reactors for Space Power Applications
The high working fluid temperature plus the gaseous nature of the fuel makes the system quite attractive for space power applications. First of all, a high working fluid temperature implies high efficiency, which in turn implies a compact system for a given power
output. A second advantage is the reduction in the size of the radiator. For any space power source, nuclear or nonnuclear, the problem of efficiently rejecting waste heat is always a major concern. Radiation is the only way of rejecting heat in space, and, for large power systems, the mass of the radiator usually tends to dominate the total system weight (Ref.2). Thus, there is a strong incentive for reducing the radiator mass. Since the rate of radiative heat transfer varies as the fourth power of temperature, the very high working fluid temperatures possible with gas core reactors offer significant advantages with regard to radiator size and, hence, system mass.
A third advantage is the elimination of the accidental criticality concerns. Before the launching of any reactor into space, extreme precautions have to be taken to make sure that the system remains subcritical throughout normal launching operations as well as during any foreseeable accidents during launch. By using a gaseous core system, the fuel in its storage container(s) is (are) in a highly subcritical state during launch. Achievement of criticality requires that the gas be both pressurized and contained or surrounded by an efficient moderator/ reflector. Also, since the process of putting a reactor in orbit is an extremely costly affair, one would like the reactor to remain in operation for a long period of time. The gas core reactor offers the possibility of long termn operation without building in a large amount of excess reactivity at the time of launch due to its capability for online refuelling and fission product removal.
5
Present Work Description
The concept analyzed in this dissertation is a gaseous core reactor power system that is designed to operate at two widely different power levels: a low power mode of operation and a burst or high power mode of operation. It is a highly enriched UF6 or UF4 fueled reactor capable of delivering power anywhere from a few MWe up to a few hundreds of MWe.
The advantage of using UF6 is that it is the only stable compound of uranium that can exist in a gaseous state at a relatively low temperature (sublimation point = 329.5 K at atmospheric pressure). This low melting point implies that the system temperature, including that of the external loop, can be kept at a reasonably low value, and still avoid condensation of the fuel on the inner surfaces of the pipes and ducts. The limitations of UF6 are that it is highly corrosive and it undergoes significant dissociation by 2000 K. As for UF4, significant dissociation starts only around 5000 K and UF4 is also much less corrosive than UF6. The problem in using UF4 as the fuel is that it is very difficult to start the system up because of the high melting and boiling points (1300 K and 1700 K, respectively, at atmospheric pressure) of UF4. One possibility is that the system could be started up using UF6 as the fuel and then once the fuel temperature gets high enough (above 2000 K), UF4 naturally becomes the fuel since one of
the major dissociation products of UF6 is UF4. A more detailed description of the system is given in Chapter 2.
The incentive for designing such a bimodal system primarily
came from Strategic Defense Initiative (SDI) requirements. The power requirements specified by the SDIO can be divided into three distinct ranges:
a) a station keeping mode requiring a few MWe for a long period of time (710 years),
b) an alert or enhanced surveillance mode requiring a few tens of MWe for shorter periods of time (a few hours to a few days) and finally.
c) a burst power or defense mode needing a few hundreds of
MWe or more for short periods of time (a few tens of minutes up to an hour or so).
The bimodal gas core reactor system that is described in this dissertation is expected to meet all of the above power requirements. i.e., the low power mode is capable of meeting the power requirements of both the station keeping phase as well as the alert phase and the high power mode of operation can handle the power requirements during the burst power phase.
As a part of the system design, both steady state (static) analysis and dynamic analysis have to be performed. This dissertation is mainly concerned with the dynamic analysis of the system.
7
Logically, the dynamic analysis can be divided into three distinct phases:
a) analysis of the system's transient behavior during the low power mode of operation,
b) analysis of the transition from the low to the high power mode of operation, and finally,
c) analysis of the system's transient behavior during the high power mode of operation.
Dynamic analysis of the low power mode has already been
performed at the University of Florida (Ref.3). This work, therefore, concentrates on the remaining two areas.
Specifically, the first part of the dissertation looks at the startup of the high power (or the transition to the high power) mode of operation.
One of the requirements which any spacebased power system for SDI burst power applications has to satisfy is the capability for rapid startup, i.e., within a short period of time the system should be able to go safely and in a predictable manner from a low power, stationkeeping mode to a high power mode. SDIO requirements limit this transition time to about 100 seconds, which is extremely short compared to the few hours taken by conventional, terrestrial nuclear plants to reach full power operation starting from low power or hot standby conditions.
The response of the gas core reactor that is analyzed in this dissertation is expected to be quite different from that of a conventional, solid core system. The highly compressible nature of the fuel, along with the fact that this is a circulating fuel system, implies that the fuel mass within the core is not going to be a constant during transients. The interaction between the core pressure, the core temperature, and the fuel mass within the core adds additional nonlinearity to the system equations. It is very important to know the behavior of key system parameters during the rapid transition from the low power mode to the high power mode of operation in order to ensure a safe transition.
To analyze the dynamic behavior of the gaseous core nuclear reactor during rapid transitions, a computer program is developed. The program (DYNAM) solves the complete system equations (nonlinear) numerically and predicts the behavior of important system parameters during the transition. It can be a valuable tool for studying the effects of variations in different system parameters on the system response and also for defining safe transition paths.
The second part of the dissertation is concerned with the
stability analysis of the high power mode of operation, i.e., assuming that a safe transition has already been made from a low power mode of operation to a high power mode of operation, this section looks at the power behavior following any perturbations (primarily, reactivity perturbations) imposed on the steady state, high power system. If the
system is found to be inherently unstable, the system design has to be modified to make it stable.
One way to perform the stability analysis is to linearize the
nonlinear system equations, and then use the powerful methods of linear stability analysis to determine the system stability. One of the objectives of this section is to determine whether a linear stability analysis is adequate or not for predicting the stability of the gaseous core reactor system. To this end, the results from the linear analysis can be compared against those obtained using DYNAM. If the linear analysis is found to be adequate, one can use it to determine quickly the range of values for reactor physical parameters and equilibrium power levels for which the system is stable, and to predict the bounds on allowable departures from the equilibrium states without having to resort to nonlinear analysiswhich can get very involved and time consuming for a system like this.
Dissertation Organization
Chapter 2 contains a review of the previous work done on
gaseous core reactor power systems. It includes a review of the gas core reactor work that has been going on at the University of Florida for the past 15 years or so. It also includes a detailed description of the bimodal reactor that is the topic of this dissertation.
This is followed by a chapter which describes the transition or startup analysis of the high power system. The models used in the computer code developed for analyzing the transition, as well as some
10
of the limitations of the code, are discussed. This chapter also presents some of the representative results generated by the code.
Chapter 4 contains the theory and a description of the methods used for the nonlinear and linear stability analysis of the burst power system. It also contains the results of the stability analysis.
The last chapter presents the conclusions. It also contains some suggestions for future work, and modifications that can be made to the computer program, DYNAM, written for analyzing the gas core reactor.
CHAPTER 2
PREVIOUS WORK ON GASEOUS CORE REACTORS
Research on gaseous core reactors has included many
independent, conceptual phase efforts for the past 30 years or so. During this time, a wide variety of system configurations and physical arrangements for the reactor core, energy conversion schemes, and fuel/working fluid form have been suggested, and a wide range of scientific feasibility studies, mostly analytical and a few experimental, have been conducted with encouraging results. As mentioned in Chapter 1, in the United States the gas core reactor concept was first suggested by Bell in 1955 (Ref. 1). He performed preliminary static neutronics calculations on a spherical core fueled with UF6 gas and surrounded by a moderator/reflector made of either heavy water (D2 0), beryllium or graphite. The impetus for much of the early work on gas core reactors came from the National Aeronautics and Space Administration (NASA). Primarily, NASA was interested in the possibility of applying such systems for rocket propulsion in space.
A number of concepts were evaluated by NASA and out of these, only two were selected for further study (Ref.2). These two were the opencycle coaxial flow system and the closedcycle nuclear light bulb concept. The general idea in both these concepts is to use the energy
12
from a fissioning plasma (of partially ionized uranium with average fuel temperatures in excess of 10,000 K) to radlatively heat the propellant which is then exhausted through a nozzle to generate the required thrust.
The original work on the opencycle coaxial flow reactor was
performed by Rom (Ref.4) and Ragsdale (Ref.5) under the direction of the NASA Lewis Research Center. The coaxial flow concept involves a spherical cavity with a low velocity inner stream of fissioning uranium metal plasma. This hot plasma is kept away from the cavity walls by a high velocity propellant gas (typically hydrogen) stream. The hot inner fuel plasma transfers thermal energy to the outer propellant stream by both convection and radiation heat transport. In this design;* the heat transfer occurs essentially unimpeded, but the drawback is that the propellant stream can physically mix with the fuel stream leading to a loss of nuclear fuel. A more detailed description of the system can be found in Schneider and Thom (Ref.6).
To avoid the mixing of the propellant and fuel streams, in the nuclear light bulb concept, the fissioning plasma is confined within a transparent (generally quartz) wall. The fuel is kept away from this transparent wall by injecting a buffer gas (argon, neon, etc.) tangentially into the cavity. Here, the energy transfer from the plasma to the propellant is by radiation heat transfer alone. Major work on this concept was performed at the United Aircraft Research Laboratory (UARL). Lantham (Ref. 7) performed a series of calculations
13
for the nuclear light bulb reactor rocket engine under the auspices of UARL.
The application of gas core reactors to terrestrial electric power generation also had some early investigations. In 1969, Gritton and Pinkel (Ref.8) of the Rand corporation reported one such study in which they evaluated a 4000 MWt spherical gaseous core system. The central gas core chamber was surrounded by a moderating/reflecting region and by banks of energy conversion devices. In another study. Williams and Shelton of the Georgia Institute of Technology (Ref.9) looked into the possibility of coupling gas core reactors with a Magneto Hydrodynamic (MHD) power conversion system for terrestrial power generation. They studied the possibility of coupling both the coaxial flow system and the nuclear light bulb system to an MHD generator. They concluded that gas core reactors, in conjunction with MHD generators, are capable of generating electricity in a safe and efficient manner. Up until then, the main obstacle preventing MHD generators from achieving their full potential was the low electrical conductivity of the working fluid. The study conducted at the Georgia Institute of Technology concluded that gas core reactors are capable of providing the high temperatures and the nonequilibrium ionization necessary for high electrical conductivities and, hence, for high MHD efficiencies.
Previous Gas Core Reactor Work At The University of Florida
As mentioned above, the major early gas core reactor research efforts in the United States were closely related to the NASA's space power needs and the program was essentially abandoned in the early 1970s. But, at the University of Florida (U of F), research on gas core reactors continued without much interruption during the 1970s and the 1980s. During these years many different concepts were studied and characterized at the U of F.
The work on gaseous core reactors started at the U of F in the early 1970s when Schneider and Ohanian (Ref. 10) proposed a nonsteady reactor called the Pulsed Nuclear Piston (PNP) engine. This was a pulsed system operating on a thermodynamic cycle very similar to the internal combustion engine with the fissionable gas as the working fluid.
This piston engine consists of a cylindrical core enclosed by a cylinderpiston arrangement just like in an internal combustion engine. This is surrounded by a moderating/reflecting region. A mixture of UF76 and helium is used as the fuel and the working fluid (the helium is added to enhance the thermodynamic and heat transfer characteristics of the primary working fluid). The fuel is introduced into the core, and is then compressed by the piston to form a supercritical system. Fission energy is rapidly released within the cylinder thereby increasing the temperature and pressure of the fuel
considerably. The fuel is then exhausted from the cylinder and the cycle is repeated. Neutronically, the core goes from a highly subcritical state to a critical state and then to a highly supercritical state.
The energy released by the fissioning gas can be extracted both as mechanical power (from the motion of the piston) and as heat from the circulating gas (e.g., through a heat exchangerturbine combination). Because of this double mode of energy extraction, the system efficiency can be as high as 50% (Ref. 11). An extensive neutronics and energetics analysis of this PNP concept was performed by Dugan (Ref. 12).
A pulsed system which is only slightly different from the PNP
was also studied by Diaz et al. (Ref. 11) in 1979. This is the Pulsed Gas Generator (PGG) concept. The PGG is very similar to the nuclear piston concept except for the elimination of the piston and making the core volume fixed. The fuel (which is again a mixture of UF6 and helium) is cyclically introduced into the core just like in a PNP. The pulsing effect is now achieved with the help of rotating absorber rods placed in the external moderating region. The PGG concept has the advantage of mechanical simplicity since it requires no moving piston. However, the efficiency attained by the system is lower than the PNP due to the fact that only heat energy (turbine power) is extracted from the system now. A more detailed description of the system can be found in Diaz et al. (Ref. 11).
A gaseous core system significantly different from the pulsed systems (PNP and PGG) was next in line to be analyzed at the U of F. This is the Heterogeneous Gas Core Reactor (HGCR). This novel concept in gas core reactor design was suggested in 1974 by Diaz and Dugan (Ref. 13). A heterogeneous core is formed by arranging bundles of moderator/coolant channels in a fissionable gaseous core region. Thus the HGCR has 'internal moderation' as opposed to the external moderation in the PM' and PGG. The core heterogeneity leads to significant thermal hydraulic and energetic advantages along with improved technological feasibility compared to other gas core reactor concepts. Ki Han (Ref. 14) performed static neutronic analysis and some kinetic analysis of this system.
Present System Description
One of the latest in the series of gas core reactor concepts that has been analyzed at the U of F is a bimodal nuclear power system. A schematic of the core proper is shown in Fig. 2. 1. Typical reactor dimensions and other key system parameters are given in Table 2. 1. The development of the system design and characterization are being done under the auspices of the Innovative Nuclear Space Power Institute (INSPI) at the U of F.
As shown in Fig. 2. 1. the system consists of two identical halves. Basically, each half of the system consists of a central fissioning region (called a Burst Power Gas Core Reactor or BPGCR) surrounded by a ring of mini cylinders (called Pulsed Gas Core Reactors or PGCRs)
Burst Power Gas Core Reactor (BPGCR) Chamber
MHD Disk Generator
Figre 1.Bimodal Gas Core Reactor Schematic
Figure 2 1.
Table 21. Bimodal Reactor Typical Dimensions and
Operating Conditions.
PGCR Data (Each Hall:
Number of PGCRs,
Height (in
Diameter (Wn
Thermal Power/PGCR (MWt) BPGCR Data (Each Hall):
Diameter (in
Height (in)
Thermal Power (MWt)
External Beryllium Thickness (in)
Fuel (UF6 + He) Mass (Kg)
Average Core Pressure (atmn) Average Core Temp1erature (K) Helium Mole Fraction in the Fuel U235 Enrichment in the Fuel
812
1.0
0.25
3.0
1.0
2.0 200.0 0.5
10.0
2030 1500.0
0.93
0.85
I
with a region of moderating beryllium or beryllium oxide in between. The PGCRs are also surrounded by a reflecting and moderating layer of beryllium or beryllium oxide. Cooling channels, carrying some inert gas like helium, are provided in the moderator region to keep its temperature below the melting point. Each half of the core exhausts fuel gas into a central duct and nozzle configuration which then directs the fuel gas into a disc MHD generator. The core fuel gas temperature is in the range of 20004000 K in this concept. This high temperature of the fuel, along with the significant amount of fission fragment induced ionization taking place in the MHD duct, is expected to provide adequate electrical conductivity and reasonable MHD generator efficiencies.
The PGCRs are identical to the pulsed gas generators described earlier. Each PGCR consists of a cylindrical core of roughly a meter in height and 20 to 30 centimeters in diameter. The core is embedded in a beryllium or beryllium oxide moderating/ reflecting region. Each half of the bimodal power system includes from 8 to 12 PGCRs, each of which is capable of delivering roughly 3 MWs of thermal power (Ref. 15). The low power mode of operation of the bimodal system is defined, for the purpose of this dissertation, as the case where only the PGCRs are supplying power. The central chamber is left empty or filled with an inert gas. The power from the PGCRs is expected to be sufficient to meet the power requirements of both the station keeping mode and the alert mode as mentioned in Chapter 1. Panicker (Ref.3)
20
has performed extensive static and dynamic neutronic analysis of the PGCR operation for this bimodal power system. Using Monte Carlo techniques, he has also studied the coupling effects among the PGCRs and between the PGCRs and the BPGCR.
The BPGCR is a cylindrical chamber about 2 meters in height and roughly a meter in diameter (each half). This core is surrounded by a beryllium or beryllium oxide region which acts as an external moderator and reflector. Each half of the BPGCR is capable of delivering 100s of MWt of power, and during the high power (or burst power) mode of operation of the bimodal system, the BPGCR will be supplying most of the system power. The focus of this dissertation is the dynamic analysis of the BPGCR.
It should be noted that since this system is only in a conceptual design state, variations of the presented design are also being examined. However, for the purpose of this dissertation, this base design is selected. The analysis that is performed in this dissertation should also be valid for other designs, provided the concept is not changed drastically.
Summaryof Previous Work Done in Related Areas
As mentioned earlier, a systematic study of gas core reactors was never undertaken in the United States. As a result of this, plus the fact that much of the work related to gas core reactors (spacebased systems in particular) was discontinued in the early 1970s, only a couple of studies on the startup and stability of gas core space nuclear
21
power systems are available. Although there are many such studies for conventional, terrestrial nuclear plants, they cannot be extended directly to the gas core space nuclear systems because of the different designs and quite different demands on these systems.
Gresho et al. (Ref. 16) and Birken (Ref. 17) deal with the startup analysis of the SNAP reactor systems (SNAP 2 and SNAP 8). These were small, solid core thermal reactors developed for meeting the power requirements of various space missions. They were capable of generating only a few kilowatts of electrical power. Since these systems were not meant for any defense applications, their startup requirements were much less severe than those of the bimodal gas core reactor system that is analyzed in this dissertation. For example, the SNAP 2 system took about 6 hrs. in going from a far subcritical state to a full power, steady state operation. This long startup time allowed the transition to be made in a much more controllable manner than for the bimodal gas core system. For the startup analysis of these systems, mostly analog computer techniques were used, and the only feedback effects considered were the fuel and grid plate temperature coefficients. The startup analysis of the gas core system is expected to be very different from the SNAP system analysis because of factors like the gaseous nature of the fuel, the excore circulation of the fuel gas mixture (which is also the coolant/working fluid), delayed neutron losses due to the fuel circulation, etc.
The SNAP 10 reactor, which belongs to the same class of
reactors as the SNAP2, was launched into space in 1965. The reactor was started up while in orbit, and it operated successfully in the space environment for more than a month. Actual experimental data are available on the startup and shutdown of this reactor. A more detailed description of the system can be found in Angelo and Buden (Ref.2).
In the startup analysis of the nuclear light bulb engine, the
authors (Ref. 18) did consider rapid startup situations (60, 300 and 600 sec. power ramps). Since the objectives of the study were to estimate the thermal stresses involved, the form of fuel best suited for the system, etc., the neutronic analysis of the startup was performed using only simplistic models. The analysis neglected delayed neutrons completely and no fuel circulation effects were taken into consideration. Also, they looked at only linear power ramps.
Lantham et al. (Ref. 19) studied the dynamic response of the
nuclear light bulb engine to perturbations induced while the system is operating at full power. To obtain the system response, they solved the full set of coupled thermal, fluid dynamics and neutron kinetics equations. 1(1 Han (Ref. 14) performed some stability analysis of the high power heterogeneous gas core reactor. His analysis was preliminary in nature and he treated fuel circulation effects in a very crude way. He also did not include the balanceofplant in his analysis.
CHAPTER 3
ANALYSIS OF TRANSITION TO HIGH POWER MODE
Introduction
The purpose of the work described in this chapter is to analyze the transition of the bimodal power system from a low power, stationkeeping mode of operation to a high power, burst mode. As mentioned in Chapter 1, this transition has to be made very fastin about 100 secondsand it is very important that key system parameters like core power, core temperature, core pressure, fuel mass flow rates, etc. be closely controlled during the transition to make sure that no maximum limits on any of these variables are exceeded. The information obtained as a result of this analysis can also be used in the overall system design. For example, the external beryllium (or beryllium oxide) moderator/ reflector has to be continuously cooled to remove heat deposited in it due to convective, radiative, and neutron/gamma heating from the core so as to maintain this region at the desired temperature. The heat deposition in the moderator needs to be known to establish the required coolant flow rates through the beryllium to prevent it from overheating and
melting. The temperature change of beryllium also introduces feedback reactivity which needs to be taken into account when determining the system's neutron multiplication factor.
For the purpose of analyzing the transition, a computer program was developed. This program can accept any userspecified reactivity variations as input, and predict the behavior of the core parameters of interest during the transition. The program is written in standard FORTRAN 77 and can run either on PCs or mainframes.
In practice, the userspecified input reactivity variation during the transition can be accomplished in a number of ways. For example, as in a conventional terrestrial nuclear power plant, the fuel mass inside the core can be kept fixed by maintaining constant and equal fuel mass inflow and outflow rates and then with the help of rotating absorber drums in the external moderator/reflector, the necessary reactivity variation can be brought about. Another way of varying the core reactivity is to vary the fuel mass inside the core. This is possible because of the gaseous nature of the fuel, plus the fact that this is a circulating fuel reactor. The fuel mass variation inside the core can be achieved in a number of different ways. One waywhich is described later in this chapter to demonstrate the working of the computer programis to control the fuel mass inflow rate into the core by controlling the fuel inlet pressure to the core. The fuel mass outflow rate from the core is then determined by the conditions existing within the core (temperature, pressure, etc.) and by the downstream
25
conditions. The code can be modified, with minimum effort, to handle any other plausible startup scenarios.
System Description
The schematic diagram of the analyzed reactor core is given in Fig. 2. 1. Typical dimensions of the core are given in Table 2. 1. For the purpose of this dissertation, it is assumed that the PGCR's are voided (unfuelled) during the transition to high power. A schematic diagram of the complete system, including the balanceofplant, is shown in Fig. 3. 1.
In Fig. 3. 1, the two halves of the BPGCR core are represented as a single core. The output from the core enters the MHD duct after passing through a convergingdiverging nozzle (not shown in Fig. 3. 1). The purpose of this nozzle is to accelerate the fuel exiting from the core to supersonic velocities before entering the MHD generator (high working fluid Mach numbers are necessary for reasonable MHD generator efficiencies). From the MHD duct, the fuel passes through a supersonic diffuser (where its velocity gets reduced considerably) and a radiator heat exchanger before getting pumped back into the core. A fuel reservoir is located before the compressor to provide any excess fuel needed while starting up the system, and to act like a sink to retain any excess fuel from the loop during burst power operational transients or shutdown. The heat exchanger transfers heat to a radiator for rejecting it to space.
Figure 3 1.
Schematic of Complete BPGCR Power System Including BalanceofPlant.
27
Referring to Fig. 2. 1, it can be seen that the reactor consists of two identical halves, and, since the PGCRs are assumed to be shut down during the transition to high power, the reactor can be essentially treated as two large identical cylindrical BPGCR's. The analysis is simplified by analyzing the behavior of only onehalf of the system during the transition. The justification for this simplification is the assumption that the neutronic coupling between the two cores is expected to be small (this coupling actually depends largely on the separation distance between the two halves, i.e., on the thickness of the separating beryllium or beryllium oxide region). The core power and mass flow rates obtained from this analysis can then be multiplied by a factor of two to get the corresponding values for the entire system.
Selection of External Moderator Thickness
A 'cavity reactor' configuration is one where the gaseous fuel
region is surrounded by a moderating/reflecting region and, unlike in a conventional solid core reactor, almost all of the neutron moderation takes place in this region external to the fuel region. The BPGCR is such a cavity reactor with a beryllium (or beryllium oxide) moderator surrounding the core. An adequate thickness of the moderating layer is crucial for achieving criticality. If this thickness is too small, the fission neutrons created within the core do not have a good chance of being moderated in the beryllium region and reflected back into the core to cause further fissions. Most of them will leak out of the system,
and, if the thickness is small enough, adequate neutron thermalization also will not occur. So, to maintain criticality, either the fuel mass inside the core has to be increased, or the core size has to be increased. At the other extreme, the thickness of the beryllium region cannot be made arbitrarily high due to system weight considerations. This is especially true in space power applications where one of the main goals of the system design is to keep the system weight down. So it is important to find an optimum thickness for the beryllium moderator.
To find the optimum beryllium thickness, the cylindrical core is first mocked up in spherical geometry by conserving the total core volume (see the next section for a discussion on the validity of this procedure). Then four group, onedimensional Sn transport theory keff calculations using XSDRNPM (Ref.20) are performed for various beryllium thicknesses. The obtained results are given in Fig. 3.2. These calculations are performed for an assumed core pressure of 50 atmospheres and a core temperature of 3000 K. The fuel atom densities can be estimated by assuming ideal gas behavior for the UF6He mixture. The calculated densities are given in Table 3. 1. As can be seen from Fig. 3.2, after the beryllium thickness exceeds about 50 cm., no significant gain in system kef is obtained by increasing the moderator thickness further. So, for the purpose of this analysis, a 50 cm. beryllium thickness is selected.
Tabe 31.Core Fuel (UF6He) Atom Densities
Nuclide Atom Density (#/cc)
U235 7.276 E+18 U238 1.283 E+18 F 5.135 E+ 19 He 1. 137 E+20
Helium Mole Fraction in the Fuel Uranium Enrichment Average Core Pressure Average Core Temperature
= 0.93 = 85 atom% = 50atm = 3000 K
Table 31.
1.20
1.00
10
Outer Beryllium Thickness (cm)
Figure 32. BPGCR keff Vs. Beryllium Moderator/Reflector Thickness
(R esut from 4 Group XSDRNPM Cations (Spherical Geometry))
Average Core Pressure = 50 atm
Average Core Temperature = 3000 K
Uranium Enrichment = 85%
0.80
0.60
0.40
0.20
20
120
Spherical Versus Cylindrical Geometr
Most of the kef calculations in this dissertation were performed by mocking up the actual cylindrical geometry with an 'equivalent' spherical geometry, and then performing onedimensional, multigroup discrete ordinates calculation using XSDRNPM. The 'equivalent' spherical geometry is formed by conserving the total core volume and the thickness of the external beryllium moderator/reflector. The reason for using a spherical geometry instead of the actual cylindrical geometry is that it is very time consuming and, hence, quite expensive to perform a 2D or 3D cylindrical geometry kef calculation on an externallymoderated gas core reactor using any standard neutronics code. One way of getting around this problem is by performing a 1 D (radial dependence only) cylindrical or spherical geometry calculation. When such 1D calculations are performed with cylindrical geometries, a transverse buckling correction factor is generally employed to account for the perpendicular (Zdirection) leakage. This buckling correction procedurewhich is routinely used in solid core reactor calculations is of no value with gaseous core reactors because of the fact that the flux within the core of a gaseous core reactor is generally very flat. A spherical geometry mockup provides a better onedimensional approximation for the cylindrical gas core reactors since the problem of perpendicular leakage is absent for this configuration.
32
To check the effectiveness of the equivalent spherical geometry mockup in predicating the kff of the actual cylindrical core, a series of calculations were performed using XSDRNPM and MCNP (Ref.2 1). The primary tool used for keff calculations in this dissertation is XSDRNPM. It is a one dimensional, discrete ordinates transport theory code from Oak Ridge National Laboratory which can perform multigroup kff calculations. MCNP is a Monte Carlo code which can perform continuous and discrete energy neutron, photon, and coupled neutron/photon transport through arbitrary 3D geometries. The MCNP calculations are taken as the standard, and the XSDRNPM results are compared against this. The XSDRNPM calculations were performed using the IBM 3090400 mainframe computer at the University of Florida and the MCNP calculations were performed using the CRAY XMP/48 Supercomputer at the San Diego Supercomputing Center.
First, a set of keff calculations were performed using XSDRNPM and MCNP on a simple, two region spherical system. The system consists of a spherical core of about 72 cm. radius and a beryllium reflector thickness of 50 cm. For this set of calculations, both XSDRNPM and MCNP treat the same spherical geometry. The results of these calculation are given in Table 3.2. An identical geometry comparison was undertaken to eliminate any geometry effects, and to pinpoint sources of error arising from neutron cross section library
Table 32. MCNP keff Versus 123Group, S4P3 XSDRNPM keff
Geometry Calculation Type keff % Difference
Spherical Continuous1
Energy MCNP 1.1140
1.01%
Spherical 123 Group 2 XSDRNPM 1.1027
1. Estimated Relative Error for MCNP Result = 0.7%
2. Flux and k~ff Convergence Criterion for
XSDRNPM Calculation = 1.0 E4 differences, including energy structure differences (see discussion below).
Calculations performed at the University of Florida, using
XSDRNPM, use a neutron cross section library which comprises data taken primarily from ENDF/BIII and ENDF/BIV compilations whereas MCNP uses mostly ENT)F/BV data maintained at the San Diego Supercomputing Center.
The MCNP code is capable of performing continuous energy neutron transport whereas XSDRNPM can perform only multigroup neutron transport. With the cross section library maintained at the U of F, the finest group structure available is 123 neutron groups. The 123group energy structure is quite accurate, and is almost as good as
a continuous energy case. So, the results given in Table 3.2 indicate that the error due to cross section library difference is quite small since the 1% or so difference shown in this table is within the statistical fluctuations of the MCNP results.
The next set of calculations were performed to get an estimate of the error involved in using a 4group energy structure compared to a finer structure, for e.g., a 26group or a 123group. Table 3.3 shows the difference in keff values obtained from a 123group, a 26group, and a 4group XSDRNPM calculations. The table also gives the CPU time taken for these three cases on an IBM 3090 mainframe. The geometry and the cross section library used are the same for all the three cases. As can be seen from the table, the 4group calculations predicts a k~ff value which is 3. 1% higher than a 26group calculation and the 26group calculation predicts a keff which is only 0.04% higher than a 123group calculation. Thus, it can be seen that the error involved in using a 4group energy structure is small enough, and the computer time savings involved is significant enough to justify the use of a 4group energy structure.
Finally, 'equivalent' spherical geometry XSDRNPM keff
calculations (both 123group and 4group) were compared against a cylindrical geometry continuous energy MCNP calculation. The 'equivalent' spherical geometry was formed by conserving the total core volume and keeping the external beryllium thicknesses the same.
Tabe 33.XSDRNPM S4P3 keff Calculations  4Group Versus
26Group Versus 123Group.
Geometry Calculation '1ype keff % Difference CPU Tme (see)
Spherical 123 Group 1.1027 311.48
0.04%
Spherical 26 Group 1.1031 21.01
3.10%
Spherical 4 Group 1.1370 0.53
1. CPU Time on an IBM
Levels of 1.0 E4.
Table 34.
3090400 Mainframe for convergence
Spherical Geometry XSDRNPM S4P3 keff Versus Cylindrical Geometry MCNP kefr.
Geometry Calculation Type keff % Difference
Spherical 123 Group
XSDRNPM 1.1027
1.2%
Spherical Continuous
Energy MCNP 1.0893
4.4%
Cylindrical 4 Group
XSDRNPM 1.1370
Estimated Relative Error for MCNP Result = 0.7%
Flux and kff Convergence Criterion for
XSDRNPM Calculation = 1.0 E4.
Table 33.
The cylindrical core is 2 meters in height and a meter in diameter with a 50 cm. of beryllium surrounding it. The 'equivalent' spherical system has a radius of 72.1 cm. and an external beryllium thickness of 50 cm. The results of the kff calculations are given in Table 3.4. There is a 4.4% difference between the MCNP result and a 4group XSDRNPM result. From the previous two comparisons (Tables 3.2 and
3.3), we know that for the same geometry, a 4group calculation predicts a keff value which is about 3% higher than a 123group calculation. In addition, we also know that the cross section library differences are responsible for small differences between the keff values predicted by XSDRNPM and those predicted by MCNP for the same geometry. So, it can be concluded that the error due to geometry difference (cylinder vs sphere) in this particular case is around 2%. This error should be smallest for "square" cylinders (i.e., cylinders with their height equal to their diameter), and should increase as cylinders become elongated.
Estimation of Gas Core Reactor Reactivity Coefficients
Using the Equivalent Spherical Geometr
Reactivity coefficients are very important neutronics parameters of a nuclear reactor, and, when designing any new type of reactor, these coefficientswhich are specific to that systemnneed to be determined carefully. They are used in the stability analysis of the system, and the reactor can be stable or unstable at certain power levels depending on the magnitude and sign of these coefficients.
37
Values for these coefficients are also needed for the startup analysis to determine the total system reactivity during startup.
Three major reactivity coefficients associated with the BPGCR are the fuel density coefficient (or equivalently, the pressure coefficient) of reactivity, the moderator temperature, and the fuel temperature coefficients of reactivity. Using the selected equivalent spherical geometry for the gas core reactor, values for these three feedback reactivity coefficients were estimated. The results are given in Tables 3.5 through 3.8.
The fuel density coefficient is the most important reactivity coefficient for gas core reactorswhere the fuel is in the form of a compressible fluid. In conventional reactors (PWRs or BWRs), where the fuel Is in solid form, this reactivity coefficient is of much less importance because of its much smaller magnitude. The estimated values of the fuel density coefficient of reactivity for the gas core reactor are given in Table 3.5. It is a prompt and positive coefficient. But it should be noted that, as a function of power or temperature, this positive coefficient translates into a negative effect, i.e., if the core power or fuel temperature increases, the density of the fuel within the core is likely to go down (since the system is an open flow system, an increase in temperature need not necessarily imply a decrease in the core fuel density) and, hence, the system reactivity will also go down.
The second reactivity coefficient associated with the fuel inside the core is the fuel (Doppler) temperature coefficient of reactivity.
Table 35. Fuel (UF6He) Density Coefficient of Reactivity
Fuel Density 1 Fractional Change in keff (gmn/cc) keff per Unit Density Change (Ak/k per gm/cc)
1.0 0.3752
2.0 0.5934+405
+ 169.50
4.0 0.8356
+ 72.79
6.0 0.9668
+ 40.78
8.0 1.0490
+ 25.50
10.0 1.1039
1. 4Group XSDRNPM S4P3 Calculation on a Spherical
Core of Radius 72.1 cm. with a Beryllium Thickness of 50 cm.
Table 3.6 shows the system kff as a function of the fuel temperature. The fractional change in keff per degree change in the fuel temperature is also given in this table. This fuel temperature coefficient of reactivity is a prompt coefficient (i.e., it follows the changes in the fuel temperature almost instantaneously) and, hence, is a very important coefficient in conventional solid core reactors from a safety standpoint.
Table 36. Fuel (UF6 He Mixture) Temperature Coefficient
of Reactivity
Fuel 1 Fractional Change in keff
Temperature keff per Degree K
(K) (Ak/k per degree K)
500 1.09867
1.766 E6
1000 1.09770
+4.191 E7
1500 1.09793
6.558 E7
2000 1.09757
1. 123Group XSDRNPM S4P3 Calculation on a Spherical Core of Radius 72.1 cm. with a Be Moderator Thickness of 50 cm.
keff and Flux Convergence Level=1.0 E5.
Increasing the fuel temperature affects the system reactivity due to changes in the fuel's microscopic cross section; this is a Doppler
effect in the fuel where resonance interaction increases due to the broadening of the resonances. In conventional reactors, the fuel enrichment is very low (34% in U235) compared to the 8085% enrichment used in the gas core reactors described in this dissertation. With such low enrichment, the fuel temperature coefficient of reactivity is dictated mainly by the Doppler effects in the huge capture resonances of the U238 nuclide. So, at all temperatures. the fuel temperature coefficient is negative. In the case of the gas core reactor considered in this dissertation (where the fuel is highly enriched in U235), the fuel temperature coefficient is dictated by very complex phenomena, where the Doppler effect in the capture resonances of U238 and U235 compete with that in the fission resonances of U235. Apart from this, there is a neutron flux spectral effect alsowhere the increased resonance interactions tend to change the neutron flux spectral distribution within the core. Some efforts were made to pinpoint the source of the oscillatory behavior of the fuel temperature coefficient of reactivity (see Table 3.6). But since no definite conclusion could be drawn from the limited study undertaken, this task is left as a suggestion for future work.
A change in the moderator temperature affects the system reactivity through two separate mechanisms. The first one is the 'spectral effect'. Since almost all of the neutron thermalization is taking place in the beryllium moderator, an increase in the moderator temperature results in a harder thermal neutron spectrum within the
41
core (due to the increased vibrational energy of the moderator atoms). A harder thermal neutron spectrum in turn implies a lower overall thermal neutron absorption (fission plus radiative capture) probability. This part of the beryllium moderator temperature coefficient of reactivity is given in Table 3.7. The values given in this table were obtained by repeating keff calculations in which only the beryllium thermal cross section data were changed, i.e., basically using beryllium cross sections with the thermal scattering kernels processed at different temperatures. Changing the thermal scattering kernel of the beryllium will alter the thermal neutron spectrum within the core and, hence, will change the neutron multiplication factor. All other parameters, including the beryllium density, were kept constant for these calculations.
As can be seen from Table 3.7, the spectral component of the beryllium temperature coefficient of reactivity switches from a small positive value at low temperatures to an even smaller, but negative, value above 600 K, and then becomes larger and more negative with increasing moderator temperature. From a safety point of view, it is desirable to have this coefficient negative throughout the entire temperature range of interest. The positive coefficient at low temperature is not much of a concern because of the fact that the temperature range where the coefficient is positive is quite small, and at normal operating power levels of the system, the beryllium temperature is expected to be significantly higher than 600 K. Also, in
42
practice, due to the large thermal inertia associated with the beryllium (or beryllium oxide) moderator/ reflector, this reactivity coefficient has a very long time constant associated with it and, hence, will contribute a negligible amount of feedback reactivity during any fast transient.
Table 37. Spectral Component of Beryllium Moderator
Temperature Coefficient of Reactivity.
Beryllium 1Fractional Change in keff
Physical keff' per Degree K
Temperature (K) (Ak/k per degree K)
296 1.08638
+3.485 E5
600 1.09789
2.186 E6
800 1.09741
1.549 E5
1000 1.09401
1. 123 Group XSDRNPM S41'3 Calculation on a Spherical Core of Radius 72.1 cm. with a Beryllium Moderator Thickness of
50 cm.
keff and Flux Convergence Level = 1.0 E5.
The second mechanism through which the moderator
temperature affects the system reactivity is the density effect, i.e., the decrease in the physical density of the moderator with increasing temperature. The beryllium moderator density coefficient of reactivity is given in Table 3.8. This coefficient has a negative effect with respect
to temperature (i.e., as the physical temperature of the moderator increases, the beryllium density will go down and, hence, the keff will also decrease).
Table 38. Beryllium Moderator Density Coefficient of
Reactivity.
Beryllium Beryllium Fractional Change Physical Density keff 1 in keff per Degree K
Temperature (K) (gm/cc) (Ak/k per degree K)
296 1.8477 1.08994
7.501 E6
600 1.8249 1.08743
9.196 E6
800 1.8063 1.08543
1.230 E5
1000 1.7863 1.08276
1. 123Group XSDRNPM S4P3 Calculation on a Spherical Core
of Radius 72.1 cm. with a Beryllium Moderator Thickness of
50 cm.
Flux and keff Convergence Level = 1.0 E5.
For obtaining the moderator and the fuel temperature
coefficients of reactivity, 123group XSDRNPM keff calculations were performed. For the fuel density coefficient of reactivity only 4group calculations were done. The reason for selecting a finer (and, hence, a more accurate) 123group calculation for the moderator and the fuel temperature coefficient determination is that the magnitudes of these
coefficients are quite small (compared to the fuel density coefficient of reactivity), and if a crude energy structure were to be used for calculating them, it may not be able to pick up the energy sensitive effects accurately enough. Also, for the same reason, very stringent flux and kef convergence criteria (105 or less) were imposed for the calculation of these two coefficients. Compared to the moderator and fuel temperature coefficients, the density coefficient of reactivity is quite large (compare Table 3.5 with Tables 3.6, 3.7, and 3.8) and, hence, only 4group calculations were performed to obtain this coefficient.
For the transition to burst power analysis, only the fuel density coefficient of reactivity is taken into account. The fuel temperature coefficient is neglected because of its small magnitude compared to the fuel density coefficient, and the beryllium temperature coefficient is neglected during the transition because of the long time it takes for the beryllium to undergo a significant temperature change. The change in beryllium temperature during the transition to burst power was calculated for various scenarios and is given later in this chapter. It was found that, even for the case where no heat is removed from the beryllium moderator, the moderator temperature change during the 100 sec. transition is not significant enough to warrant the inclusion of the beryllium temperature coefficient of reactivity in the transition. In contrast, when the system is operating at high or burst power levels, and full power transients are analyzed, the moderator temperature
coefficient is one of the significantalthough slowfeedback effects. So, in Chapter 4, where the stability analysis of the system is described, this coefficient is included.
Models Employed in the Computer Program Modelling of the Core Geometr
As mentioned earlier, the PGCRs are assumed to be shutdown
during the transition to high power. This makes it possible to treat the reactor as two identical halves (BPGCRs). The analysis is further simplified by assuming a weak neutronic coupling between the two cores and, hence, only one half of the system is treated during the transition. One half of the BPGCR is shown in Fig. 3.3. The core is taken to be 2 meters in height and a meter in diameter. It is surrounded by a 50 cm. thick beryllium moderator/reflector (not shown in the figure).
As shown in Fig. 3.3, the BPGCR is modelled as an 'open flow' system with no valves at the core inlet or exit; in a real system, there may be some sort of valving or orificing at the core inlet for controlling the fuel flow into the system and, if necessary, for shutting the system down. Similarly, there might be a valve or orifice at the core exit also to control the fuel flow out of the core. The flow characteristics will be altered drastically if the flow through the core is controlled in this manner and, hence, the fluid flow analysis performed here would have to be modified to take this into account. A converging nozzle is provided at the core inlet. At the core exit there is a
BPGCR
Figure 33 Single Cavity BPGCR
Wnet
Exit '
convergingdiverging (CD) nozzle. The converging nozzle at the core inlet is expected to provide an 'inherent' control mechanism while the system is operating at full power. For example, if the system power level goes up for some reason, the core pressure and temperature will also go up, and this will cause a reduction in the fuel inflow rate (provided the flow at the inlet is nonchoked and the core inlet pressure is kept constant). The reverse action occurs if the system power drops for some reason. The CD nozzle at the core exit also will provide a similar controlling action. The exit nozzle is assumed to be operating under choked conditions. The CD nozzle is also needed at the exit for accelerating the fuel to supersonic velocities for the MHD generator. The arrangement of the balanceofplant is shown in Fig.
3.1.
Thermodvnamic and Heat Transfer Modelling
The fuel (UF6 He mixture) is assumed to behave like an ideal gas with constant specific heats (for the purpose of this study, changes in the chemical composition of the fuel are neglected). A lumped parameter approach is selected for the thermodynamic and heat transfer treatment of the system. For example, it is assumed that an average value of the fuel temperature may be used to describe the heat transfer properties of the fuel inside the core at any given time. The same is true for thermodynamic properties like pressure, enthalpy, etc.. An average value is also assumed for the temperature and the heat capacity of the moderator.
48
So, for example, we can write differential equations governing the dynamic behavior of the average core temperature from a simple heat balance of the form
Core Heat Accumulation Rate = Core Heat Deposition Rate
 Core Heat Removal Rate. (31) The rate of heat accumulation in the fuel is given by
where
Mf Mt = instantaneous fuel (UF6+He) mass inside the core (Kg),
Cpf = average fuel specific heat at constant pressure (JfK gK),
and
Tf Mt = average instantaneous core fuel temperature (K).
The rate of heat deposition within the core is taken to be a constant fraction of the instantaneous core fission power, i.e.,
Qcore = 0.9 * N (32) where
Qcore = heat deposition rate within the core (watts), and
N = instantaneous core fission power (watts).
So, it is assumed that ninety percent of the total fission energy generated within the core gets deposited within the core. The other ten percent escapes from the core in the form of prompt and delayed gammas plus the kinetic energy of the fission neutrons. Almost all of this energy is assumed to be deposited in the moderator/ reflector.
Three mechanisms of heat loss from the core are considered:
radiative transfer to the wall/moderator, convective heat transfer, and the loss due to the energy transported by the flowing fuel into and out of the core.
The net radiant heat transfer rate between the core and the
wall/moderatorwhich are taken to be at uniform temperaturecan be written as (see Ref. 22 for details)
qrad = Ac G (e Tf4  a TW4) (33) where
qrad = net radiant heat transfer rate between the core
and the wall (watts),
A= surface area of the wall/beryllium moderator which is
exposed to the core (fuel) (in2),
G StefanBoltzmann's constant (watts / m2K4) Tf =average fuel temperature inside the core (K),
Tw= average wall/moderator temperature (K),
P_ = average fuel emissivity, and
a = average fuel absorptivity.
(Implicit in the use of the above equation is the assumption that the two bodies are at uniform temperature, i.e., radiation from each body is characterized by a single average temperature.)
50
If a further simplifying assumption of a = F, (i.e., Kirchhoff's law) is made, the above equation becomes
qrad =A (TF,(Tf'  T.(34 ) (Note that in actual practice, radiative heat transfer from a gas presents a far more difficult problem than is indicated by the above equation. It is complicated by the spectral and directional effects involved with gaseous absorption and emission, and also by the dependence of both (X and E on temperature. But, for the present work, the simplified treatment given above is presumed to be sufficient.)
The film coefficient of heat transfer, H, for incompressible pipe flow is generally based on the difference between the wall temperature and mean stream temperature. For compressible flow at high speeds, it is more appropriate to base H on the difference between the actual wall temperature Tw, and the 'adiabatic' wall temperature, Taw (see Ref.23 for details). Accordingly, the total convective heat transfer rate between the core and the wall/reflector can be written as
qconv =H A(Tw Taw) (35) where
qconv = net rate of convective heat transfer (watts),
H =effective heat transfer coefficient from fuel to the wall
(watts/rn2 K), and
A =surface area of the wall/beryllium moderator which is
exposed to the flowing fuel in the core (in2). The adiabatic wall temperature is defined by
Taw = 1 + R  1 M
T 2
where
T = average static temperature of the stream (K),
R = recovery factor (value between 0.87 and 0.91 for subsonic
pipe flow),
M =Mach number of the flow, and
y =ratio of specific heats for the flowing gas (fuel).
The film coefficient of heat transfer, H, is defined for turbulent flow (the flow through the core will be highly turbulent even for relatively small fuel velocities of the order of 20 m/sec) with the help of the dimensionless Nusselt number, Reynolds number and Prandtl number as follows:
Nu = 0.0364 (Re * PrjO.75 (36) where
Nu =Nusselt number,
Re =Reynolds number, and
Pr =Prandfl number.
(See Ref.23 for details of the above treatment of convective heat transfer.)
Neglecting the changes in the kinetic and potential energies of the flowing fluid across the core, the energy loss due to the energy transported by the flowing fluid into and out of the core can be written as
qtrans = mout Cpf Tout  mi Cpf Tin (37) where
mout (t) = instantaneous fuel mass flow rate out of the core
(Kg/ sec),
MO= instantaneous fuel mass flow rate into the core
(Kg/ sec),
Tn= Average fuel temperature at the core inlet (K, and
Tot= Average fuel temperature at the core exit (K).
This includes the energy loss associated with the thermal energy changes across the core plus the flow work changes.
Equations 3.2, 3.4, 3.5, and 3.7 can be substituted in Eq. 3.1 to get a final expression for the average core fuel temperature.
A similar treatment is also performed for finding the
instantaneous average beryllium moderator temperature, i.e.,
Be Heat accumulation rate = Be Heat deposition rate
Be Heat removal rate (38)
The rate of heat accumulation in the beryllium moderator is given by
MBe CP Be (dTBe / dt)
where
MBe = mass of beryllium moderator/ reflector (Kg),
CpBe = average specific heat of beryllium (J/KgK), and
Te= average beryllium temperature (K).
The rate of heat deposition in beryllium has three components to it. The first one is the heat deposition due to the neutrons and gammas escaping from the core. This is taken to be a constant fraction of the instantaneous core fission power. So,
QI~e = 0.1 * N (39)
where
QBe = rate of heat deposition due to neutrons and gammas in
the beryllium (watts), and
N = instantaneous core fission power (watts).
(That is, here it is assumed that the ten percent energy lost from the core due to escaping neutrons and gammas is given up fully in the moderator. See also Equation 3.2.)
The second source of heat deposition in the beryllium
moderator is the radiative transfer from the core. This part is given by Equation 3.4 above. The third source of heat transfer to the beryllium is the convective transfer from the flowing fuel inside the core to the beryllium wall. This contribution is given by Equation 3.5.
External Loop Modelling
An ideal Brayton cycle is taken to be the thermodynamic cycle. The TS diagram for the cycle is shown in Fig. 3.4. The heat addition process within the core, and the heat rejection to the heat exchanger/ radiator are assumed to take place at constant pressure. The MHD generator is modelled as a conventional turbine. The performance of the turbine is then dictated by selecting an enthalpy extraction ratio and an isentropic efficiency for the turbine (these are specified as constants to the program). The values selected for the enthalpy extraction ratio and the isentropic efficiency are based on some earlier analysis of the MHD generator performed by Gerard Welch at the U of F. The enthalpy extraction ratio, together with the isentropic efficiency, can be used to obtain values for the stagnation pressure and temperature change across the turbine (MHD generator). The performance of the compressor is handled by specifying an isentropic efficiencywhich is also taken to be a constant.
In general, the analysis of the heat exchanger/radiator is very involved and, in the absence of detailed prior analysis, it cannot be treated in a simplistic way as is done for the MHD generator or the compressor (see discussion above). That is, the behavior of the heat exchanger/radiator cannot be described by one or two performance indices as is done for the compressor or the MHD generator.
To place the treatment of the heat exchanger/radiator at the
same level as the analyses of the MHD generator and the compressor,
(stagnation
values)
I
p =p 02 03
TURBNE
0O4 P0 1
Wout
Figure 34. TS Diagram for the BPGCR Power System
56
and also to reduce the computational time for the transition to burst power analysis, one of two alternate simplifying assumptions can be made: either it can be assumed that the heat exchanger/radiator combination rejects heat at a constant rate or it can be assumed that it rejects a constant fraction of the heat generated within the system (the total heat includes the heat generated within the core due to fissions occurring there plus the work done by the compressor on the system). For a gas thermodynamic cycle like the Brayton cycle, the work done by the compressor is quite significant (see discussion below). Note that selection of either assumption implies that the radiator has to be of variable size (area). This is because of the fact that during the transition to burst power, the heat generated by the system (and, hence, the temperature of the working fluid) is continuously increasing, and if it assumed that the heat rejection rate is constant or is a constant fraction of the heat generated within the system, then the radiator size must be changing continuously. Such variable area radiator concepts have been proposed.
For the purpose of this analysis, the second alternative is
selected since it predicts an increasing heat rejection rate during the transition and, hence, is closer to the actual situation than the first assumption (i.e., in reality as the system power increases, the amount of heat rejected is also expected to increase, and dramatically so because of the T4 dependence of the radiative heat transfer).
Brayon Cycle vs Rankine Cycle
For the purpose of this dissertation, the thermodynamic cycle is taken to be an ideal Brayton cycle. The primary motivation for employing a gaseous cycle like the Brayton cycle was that the fuel/working fluid considered initially was a UF6 He mixture, and at typical operating temperatures and pressures existing within the system, this mixture always remained in a gaseous form. Also, initially, the possibility of a pulsed system with acoustic amplification was considered for the burst power mode of operation. It was conjectured that under these conditions, a Brayton cycle might prove to be satisfactory. INSPI has also looked into the possibility of using a fuel/working fluid mixture consisting of UF4 and helium. The UF4 can be condensed out of the mixture at the temperatures and pressures existing within the system and, hence, the possibility of a hybrid Rankine/Brayton thermodynamic cycle was also considered. Currently, INSPI is investigating fuel/working fluid mixtures of UF4 and metal fluorides in which the entire fuel/working fluid mixture is condensible under expected operating conditions. Such a fuel/working fluid mixture means that the system can operate on a true Rankine cycle.
The main drawback of the Brayton cycle is the compressor work needed to compress the gas. This is significantly higher for the Brayton cycle than for a Rankine cycle. In the latter, since only an almost incompressible liquid is pumped, the typical backwork ratio
58
(which is defined as the ratio of the pump or compressor work to the turbine work) is only a few percent or less compared to 5080% for a Brayton cycle (for an ideal Brayton cycle it is of the order of 50% but in actual systems, it can be as high as 80%).
A second drawback of the Brayton cycle, which is very significant for space power systems, is the temperature at which heat rejection occurs. In a Rankine cycle, the heat rejection can take place at constant temperature (condensing vapor) whereas in a Brayton cycle it occurs at a varying temperature. If we were to determine an "effective" constant temperature for heat rejection for a Brayton cycle, it turns out to be much lower than the constant temperature at which the Rankine cycle rejects heat (provided we place comparable thermal limits on the Brayton and Rankine cycles). This factor can dramatically increase the size of the radiator.
If the frictional losses in the system are neglected, and if the efficiencies of a Brayton cycle and a Rankine cycle are compared (for roughly the same operating conditions), the Rankine cycle can actually turn out to be poorer than a Brayton cycle (2226% for a Brayton cycle compared to roughly 22% for a Rankine cycle) (numbers are based on some earlier work done by Gerard Welch at the U of f)). However, we cannot neglect frictional loses, and making up for these in a Rankine cycle costs almost nothing; the cost in a Brayton cycle can be very high because of the high backwork ratio. With the frictional loses included, the Brayton cycle efficiency turns out to be lower than that for the
Rankine cycle. This lower efficiency also translates into a bigger radiator mass (since a lower efficiency implies more heat rejection for the same electrical power output). Fluid Flow Modelling
As mentioned earlier, the fuel (UF6 He mixture) is assumed to behave like an ideal gas with constant specific heats. The fuel is also assumed to maintain its integrity throughout the transition to burst power (i.e., no dissociation of the UF6 is taken into account).
Even after the commencement of the transition, it takes some time before the core temperature and pressure builds up to levels where significant dissociation of UF6 can occur. Significant UF6 dissociation occurs only at temperatures above 1600 K at atmospheric pressure. With increasing gas pressure, dissociation becomes significant only at higher temperatures. For e.g., with a UF6 partial pressure of about 10 atm., the fuel can exist without much UF6 dissociation up to about 2200 K. The fuel temperature profile during the transition is given in Fig. 3.21 later in this Chapter. From this figure, it can be seen that, at least for the first 50 sec. of the total 100 sec. transition time, the fuel temperature remains below the 2200 K value. (UF6 may be able to exist without significant dissociation at much higher temperaturesup to about 3000 K, provided fluorinating agents like ClF3 or BiF3 are added to the fuel (Ref.24)).
If the fuel is in the form of UF4, the dissociation is not at all a
concern even at temperatures as high as 30004000 K. Neutronically, it doesn't matter whether the fuel is in the form of UF4 or UF6. The extra fluorines have no significant neutronic effect. So, the neutronic analysis done here remains valid whether UF6 or UF4 is used as the fuel.
Core hydrodynamics is represented by fully developed,
onedimensional, turbulent flow (as mentioned earlier, the fuel flow through the core will be highly turbulent even for relatively low fuel velocities of the order of 2050 m/sec). Frictional effects are neglected. Isentropic flow is assumed for the flow through the inlet and exit nozzles so that stagnation temperatures and pressures of the fuel will adequately describe the flow of the fluid through the nozzles. For example, the maximum fuel flow rate through the exit nozzle is given by the relation
hmax = 0.6847 PO A
R (T)"/2 (310) where
mmax = maximum fuel mass flow rate (Kg/sec),
A = exit nozzle throat area (in2),
Po = core stagnation pressure (Pa), and
To = core stagnation temperature (K.
(This maximum flow occurs under choked flow conditions, i.e., when the Mach number at the throat is unity). The exhaust nozzle is assumed to be under choked flow conditions always (flow through the convergingdiverging nozzle at the core exit has to be choked to obtain supersonic flow going into the MHID generator). The converging nozzle at the core inlet is operated under nonchoked conditions for reasons given earlier in the geometry modelling section. Neutron Kinetics Modelling
It was mentioned earlier that it is assumed that negligible
decomposition of UF6 takes place during the transition to burst power. Apart from this, it is also assumed that the fuel composition remains constant while going through the external loop. This implies that the UF6 mole fraction of about 93% and the 85% uranium enrichment remains the same. It is also assumed that there is no significant buildup of strongly neutron absorbing fission products. These assumptions should not result in too much error since the fuel burnup and fission product buildup are going to be very small during the short transition time.
To model exactly the time behavior of the neutron population within a reactor, use of the timedependent Boltzmann's neutron transport equation has to be made. This is all the more true for the type of gas core reactor that is analyzed in this dissertation because of the fact that within the core proper, almost no neutron moderation
(scattering) is taking place, and the core is only of the order of a mean free path or less in size. The only significant interaction taking place within the core is the absorption of thermal neutrons returning to the core from the reflector/moderator. But, from a practical standpoint, implementing this model can be very tedious and computationally intensive, and is attempted only in major time dependent transport theory computer codes.
Some previous studies conducted at the U of F indicate that, for gas core reactor systems which posses a significant degree of grayness, multigroup diffusion theory is generally adequate for static calculations provided good fast and thermal group constants are used (Ref. 25) (grayness for an externally moderated gas core reactor is defined as the ratio of the thermal neutron current into the core from the moderating/reflecting region to the value of the thermal neutron flux at the coremoderator/reflector boundary). For small transients (i.e., small reactivity variations), we can assume that timedependent, multigroup diffusion theory is adequate to describe the dynamic behavior of the gas core reactor. Sometimes, even this model is too detailed for practical implementation in reactor kinetics analysis due to excessive computational requirements.
For this reason, in this dissertation (as is the common practice for many preliminary level reactor kinetics studies), it is assumed that the neutron dynamic characteristics of the gas core reactor can be adequately described by the average neutron density, 'n', i.e., the
spatial dependence of the neutron population is neglected. Often, in neutron kinetics studies, 'n' is assumed to represent only the neutrons in the energy range where most of the fissions occur for example, thermal neutrons in a thermal reactor. This type of analysis is by no means exact but, for preliminary scoping purposes, this is adequate. It is capable, for example, of providing the time dependence of the gross reactor power level during transients.
With the above assumption, the dynamic behavior of the neutron density within the core can be represented by
dn(t) _P(t)  feff
dt 1' n(t) + Xi Ci(t) (.1 where
n(t) = core neutron density at time t (#/cm3),
p (t) = time dependent system reactivity,
Pff= effective delayed neutron fraction, V=mean neutron generation time (sec),
ct(t) = effective delayed neutron precursor density for the ith
delayed neutron group (#/CM3), and
X1=decay constant for the ith delayed neutron precursor
group (sec1).
To complete the description of the neutron dynamics of the
reactor, an equation for the ith delayed neutron precursor density, c1, also has to be provided. For the gas core reactor analyzed here, this
precursor density equation is quite different from that for a conventional, solid core reactor. Since this is a circulating fuel reactor, not all of the delayed neutron precursors will decay within the core and thereby contribute to the total core neutron population. Some of the precursors are going to be swept out of the core (along with the fuel) before they can decay. Now two things can happen to these precursors which have left the core without decaying: they can either decay outside the core (so these neutrons are lost forever, and will not contribute to the neutron population within the core) or they can reenter the core after going through the external loop without decaying. In either case, the effective fraction of delayed neutron for this type of reactor is going to be less than that for a stationary fuel reactor. The precursor density equation taking into account the fuel circulation is
dci(t) = IPeff n(t) x.1 c1(t)  c1(t) + c1(tx1l) exp (x1 X)i
dt I1* T (312)
for i = 1,2,...,6.
where
=c mean core residence time for the fuel (see), and
Tj mean external loop circulation time for the fuel (sec) (All other variables are as defined in Equation 3.11.) The third term on the right hand side of the above equation accounts for the loss of the delayed neutron precursors from the core due to
flow, and the last term on the right hand side takes care of the reentry of the precursors at the core inlet after going through the external loop.
So, the effect of the fuel flow outside the core is to make the
reactivity needed (and, hence, the fuel mass needed inside the core) for steady state operation larger than the value needed if it were a nonflow system.
In Equation 3.11, 'n' represents the instantaneous neutron
density within the core. But it can represent any integral or volume averaged property of the reactor that is proportional to the instantaneous neutron density. For example, it can be redefined so that it represents the total number of neutrons within the core at any given instant or the fission rate or the total core power at any given instant, etc.. For convenience, the dependent variable in Equation 3.11 is defined to be the instantaneous power, N, of the reactor. Now the precursor density, c1, used in Equations 3. 11 and 3.12 has to be replaced by a new quantity C1 defined as follows:
Cj = WfV YXf c1 V (313)
where
C= effective delayed neutron precursor density of the ith
delayed neutron group (#/cm3),
Wf = usable energy produced per fission (MeV/fission),
If =average thermal macroscopic fission cross section (cmi1),
V =total core volume (cm3), and
v =neutron velocity (cm/sec).
As a final simplification, instead of the usual six delayed neutron groups, only one effective group of delayed neutrons is considered, i.e., the summation term in Equation 3.11 is replaced by one term which represent an "effective" delayed neutron group. Also, Equation 3.12which is a set of six equations, is now replaced by one equation for the "effective" delayed neutron group. This simplification is also done so as to reduce the computational time. The computer program developed can equally well handle all six delayed groups (with a corresponding increase in the computational time).
There are different ways of defining an "effective" one delayed group of neutrons. The recipe used in this dissertation for obtaining the value for 6, and the corresponding value for X is as follows;
5~X~ (314) and
(XP~ X~)(315)
for i =1, ..6
where
_ = fraction of the ith delayed neutron group emitted per
fission, and
X= decay constant for the delayed neutron precursor for the
ith delayed neutron group.
The neutron generation time, 1% used in Equations 3.11 and 3.12 is an important factor, and is defined as the average time between the birth of a neutron and its subsequent absorption inducing fission. It is related to the prompt neutron lifetime, T1 of the reactor by the relation
1P = 1I /kef
where 'kf is the effective neutron multiplication factor of the reactor. The prompt neutron lifetime (which is defined as the average time between the birth of a neutron and its subsequent removal from the systemeither by absorption or by leakage) is a very difficult parameter to estimate analytically since it depends in a complex manner on a number of system variables like geometry, fuel enrichment, reflector thickness, etc. So, for the purpose of this work, T1 for the system is determined by independent Monte Carlo calculations, and the value thus obtained is used as an input to the program.
Programminif Aspects
The computer program, DYNAM, developed as a part of this
work for analyzing the transients, is written in standard FORTRAN 77. It can run either on a PC or a mainframe (the program calls for certain timing subroutines which are only available on IBM mainframes. These call statements need to be deleted from the program before trying to run it on a PC or on a nonIBM mainframe).
Simple finite difference techniques are used for solving most of the relevant equations like those for core temperature, beryllium temperature, core fuel mass, etc.. The coupled neutron kinetics equations are solved using a fourth order RungeKutta method (this method is selected because it is quite easy to implement, plus it has good numerical stability characteristics).
The only difficult term to handle in the neutron kinetics
equation is the precursor reintroduction term, i.e., the last term on the right hand side of Equation 3.12which represents the precursors which are reintroduced into the core after going through the external loop. Thbis term is complicated because, at any given instant, it depends on the precursor density which had existed within the core at a time Ti (the mean external loop residence time) before.
When the core power is varying rather slowly, a reasonable appro.x imation would to be take
This is certainly not true during the transition to burst power where the power goes from almost zero to a few hundred megawatts in about 100 sec. So, in DYNAM, the precursor reintroduction term is handled numerically as described below.
The precursor density at each time step is stored in a stack,
and, at any instant, the value of 'C' which existed at a time tj before is used to estimate the precursor reintroduction term. Once a value C(tt1) is used during a time step, it can be discarded from the stack of stored values of 'C'. That way, at any given instant, the number of values of 'C' that need to be stored is only
tri / step
where 'step' is the (constant) time step size used in the numerical computation.
Sometimes this can create memory limitation problems while running on a PC. In this case, the precursor density can be averaged over a fixed number of time steps before storing it. The number of times steps over which the precursor density is averaged depends on the transient. If it is a very fast transient, one might want to average the density over a small number of time steps, whereas, for slower transients, the density can be averaged over a large number of time steps without incurring serious error.
Testing of the Computer Program
Before the developed computer program, DYNAM. can be used for analyzing the transition to burst power, it has to be tested or benchmarked. The complete program which includes the neutron kinetics, the thermodynamic, the heat transfer, and the fluid flow aspects of the gas core reactor, cannot be benchmarked against any standard codes since there are no eisting codes which are capable of analyzing transients for such a novel type of reactor (novelty includes the combination of gaseous fuel in a closed circulating loop, highly compressible fuel which also serves as the working fluid, etc.). So, it was decided to test individual sections of DYNAM against some standard codes or against analytical results where possible.
One of the major subroutines in DYNAM is the one which
estimates the effective neutron multiplication factor of the system (or the system reactivity). This subroutine accepts the reactor dimensions (core radius and beryllium thickness) and the fuel atom densities as input, and it estimates the system kff. The system keff needs to be known at every instant during the transition to burst power.
Since the system keff has to be estimated a large number of
times during the transition analysis, this subroutine was written giving speed of calculation priority over accuracy. The subroutine mocks up the cylindrical core with a spherical geometry, and then performs a tworegion (core and moderator/ reflector), twogroup (fast and
thermal) diffusion theory keff calculation. The microscopic group constants needed for the keff calculation are provided as input to the subroutine. These group constants are obtained from independent. twogroup transport theory calculations using XSDRNPM. These microscopic constants are found to be not very sensitive to small changes in fuel atom densities within the core or to changes in the core dimensions. But, if either of the two (geometry or fuel atom densities) is changed significantly, the group constants need to be reevaluated using XSDRNPM. The effects of temperature on the microscopic constants are handled via feedback coefficients. These coefficients are given in Tables 3.6 and 3.7.
The system keff values calculated by the DYNAM subroutine are compared against those estimated by XSDRNPM (for the same geometry and fuel atom densities). Comparisons are made against twogroup diffusion theory and twogroup transport theory results from XSDRNPM. The results of the comparisons are shown in Table 3.9. The comparisons are done for various fuel atom densities within the core.
As can be seen from the table, in all cases, the error involved is less than a percent and, hence, is quite acceptableespecially considering the fact that the subroutine is quite fast. It takes less than
0.01 sec. to estimate the keff of the system to a convergence level of better than 104 on an IBM 3090 mainframe. Compared to this, a
Comparison of keff Values Obtained Using the keff Subroutine of DYNAM and the XSDRNPM Code for Different Fuel Density Factors.
1. A fuel density factor of unity corresponds to the
densities (in atoms/barncm);
U235 U238
F
He
following atom
1.458 E5 2.548 E6 1.027 E4 2.275 E4
2. Two group, two region, diffusion theory calculation on a
spherical core with 72.1 cm. core radius and 50 cm. beryllium
moderator/ reflector thickness.
3. Two group, two region diffusion theory calculation on a
spherical geometry of the same dimensions as given above.
4. Two group, two region transport theory calculation (S4 P3)
on a spherical system of the same dimensions as given above.
Table 39.
keff Percentage Difference Core 3Between DYNAM keff
Fuel 2Group 2Group and____ ____Density2 XSDRNPM XSDRNPM
Dentyr DYNAM (Diffusion (pd Two Group Two Group
FactorTheory (P XSDRNPM XSDRNPM Option) (Diff. Theory) (S4P)
0.2 0.6990 0.7025 0.7000 0.498 0.143 0.5 1.0356 1.0411 1.0384 0.528 0.269 0.8 1.1760 1.1789 1.1758 0.246 +0.017 1.0 1.2310 1.2337 1.2302 0.219 +0.065
twogroup, tworegion XSDRNPM kff calculation on a similar system takes about 0.5 sec.
The fact that the DYNAM results are closer to the two group
transport theory results is purely coincidental. In general, it was found that diffusion theory predicts a kff which is higher than the keff prediction of transport theory. So, the XSDRNPM two group transport theory results are lower than the XSDRNPM diffusion theory results. It so happens that DYNAM predicted keff values are also slightly lower than the XSDRNPM two group diffusion theory results and, thus, actually a little closer to the two group transport theory results.
To complete the benchmarking of the keff subroutine of DYNAM, the fast and thermal group fluxes obtained using the kff subroutine are also compared against those obtained using two group transport theory calculations (XSDRNPM). The results of this comparison are given in Figs. 3.5 and 3.6. As can be seen from the figures, except for the case of fast flux within the core, DYNAM predictions of the thermal and fast fluxes are in agreement with the predictions of XSDRNPM. DYNAM underpredicts the fast flux within the core by roughly 10%. The fluxes shown in these figures (group 1 and group 2) are normalized such that the total fission neutron source in the core is one neutron per second (it should be mentioned here that the agreement between flux shapes predicted by the diffusion theory and the transport theory is not always this good. For certain fuel density regimes, the agreement
2.OOE4
0
Figure 35.
20 40 60 80 100 120 140
Radial Distance From the Center (cm)
Comparison of Radial Flux (Fast) Profile Within the BPGCR
U
Figure 36.
20 40 60 80 100 120 140
Radial Distance From the Center (cm)
Comparison of Radial Flux (Thermal) Profile Within the BPGCR.
between Sn, and diffusion theory fluxes can be very poor even though the kff values are in reasonable agreement).
A second major subroutine in DYNAM is the one which solves the coupled point reactor kinetics equations. This subroutine uses a fourth order RungeKutta method to solve the coupled neutron and precursor density equations. Two separate methods were used to check the validity of this subroutine. First, the effectiveness of this subroutine in predicting the neutron density behavior after step insertions of reactivity was compared against analytical results. For two different step reactivity insertions, the results obtained using the developed subroutine are compared against the "exact" analytical results (see Ref. 26), and also against the 'prompt jump' approximation results (the prompt jump approximation assumes that immediately following a step change in reactivity, the delayed neutron precursor concentration remains unchanged over the time of the sudden rise or drop in flux). The "exact" analytical result is for one delayed neutron group; it gives good results when the reactivity step is no greater than about 0.3 3 (see Ref.26 for the limitations of this "exact" analytical expression). Figures 3.7 and 3.8 gives the results of this comparison. Again, it can be seen that the results from the developed subroutine agree quite well with analytical predictions.
Next the developed subroutine is compared against the standard point reactor kinetics code ANCON (Ref.27). For the purpose of
0 0.4 0.8 1.2 1.6 2
Time After Step Insertion of Reactivity (sec)
Figure 37.
Benchmarking of DYNAM Against Analytical Methods  Positive Reactivity Insertion Case
4A ~ 0.0006 sec
4.1
0
A
e ROP DO 00. . . .
TieAtrSe0Isrinoeatvt sc
Fiue38Ce6mrigo YA gas nltclMtosNgtv ReciiyIsr0nCs
comparison, three different ramp reactivity insertion rates were selected. Either six delayed neutron groups or one delayed neutron group treatment can be used in ANCON. Figures 3.9 and 3. 10 compare the results from DYNAM against six delayed neutron group ANCON results (note that the one delayed group ANCON results agree exactly with the DYNAM results).
Since ANCON cannot deal with circulating fuel type reactors (it assumes that the fuel is stationary), for the purpose of comparing the results of DYNAM with those of ANCON, the fuel circulation effects (for example, the loss of delayed neutron precursors from the core, the delayed neutron precursor reintroduction, etc.) are zeroed in DYNAM.
Results
The developed computer program, DYNAM, is now used to study the transition to burst power. As mentioned before, the program is very versatile and can be used in its present form or with very little modification to study different aspects of the rapid transition. In the following sections, a few selected capabilities and applications of the program are demonstrated.
First, two separate reactivity insertions are applied to the core to bring it up from almost zero power to full power in roughly 100 seconds. During this transition, crucial system parameters like core pressure, core temperature, etc., are estimated by DYNAM. Based on the behavior of these parameters, system design can be developed or modified. For example, if one finds that, for some desired operating
1011. 10 0
10 8 10 7 10 6 10 5 10 4 10 3 10 2 10 1 10 0
0 4 8 12 16 Time After Ramp Reactivity Insertion (sec)
Figure 39. Benchmarking of DYNAM Against ANCON  Case 1
10 7'
9k ANCON (6 Gp)
10 6. 0 DYNAM (1Gp)
Ii) & ANCON (I Gp) 10 5Q
104
Rho Wt = +0.001 t 4' N0)= 1.0 S 1030 1= 0.0065 4' 1P = 0.0006 sec
102.
101.
100
0 2 4 68102
Time After Ramp Reactivity Insertion (see)
Figure 310. Benchmarking of DYNAM Against ANCON  Case 2
power level, the average core temperature is higher than some acceptable value, the system design can be modified, e.g., to increase the fuel flow rate through the core so as to reduce the average core fuel temperature. This increased fuel flow rate through the core might require changing the ratings of the pumps or even changing the core inlet and eit nozzle areas. The behavior of these crucial system parameters can also be used in selecting the best possible reactivity insertion pattern that is capable of taking the system from zero power to full power in a safe manner.
The way DYNAM is set up currently, to bring about the transition to burst power, one has to specify the desired reactivity variation (this can be almost any conceivable reactivity variation). The program will then estimate the fuel mass variation needed within the core to bring about this reactivity variation. It also determines the fuel flow rates needed in order to obtain the necessary fuel mass inside the core, and the core inlet pressure is automatically adjusted to give the required fuel flow rates.
There are many types of reactivity variations that can be applied to take the core from zero power to full power. For the purpose of demonstration, two different reactivity insertions are selected. The first one is a pure linear (or ramp) reactivity variation. In fact, as shown in Fig. 3.11, it is a combination of two linear reactivity variations. To simulate the conditions in an actual transition scenario, the core is initially assumed to be far subcritical (ptot = 1 .0). This
(Figure Not to Scale)
1.0
0 20 40 60 80 100
Time After Startup Intiation (sec)
FIG. 3.11 BPGCR Input Reactivity Variation With Time  Case I
would then correspond to a situation where the fuel is just being introduced into a previously unfuielled (voided) central chamber.
During the first part of the reactivity variationwhen the core is subcritical, a very steep reactivity variation of 0.02 Ak/k per second is introduced into the core (for a U235 fueled system with a delayed neutron fraction of about 0.0065, this number translates to roughly $3.0 per second). Once the system reactivity becomes zero (keff = 1.0), the slope of the linear reactivity variation is reduced to 0.000 15 Ak/k per second (roughly $0.02 per second).
Figures 3.12 through 3.16 show the variation of some of the
important system parameters during the rapid transition to about 100 MWt as computed by DYNAM. The core power is set at an initial value of 1000 watts within the program. This might correspond to a situation where the bimodal system is getting ready to go to a high power mode of operation from a low power mode, and some power is generated within the BPGCR due to neutrons reaching it from the surrounding PGCRs. The behavior of the core power (Fig. 3.12) can be explained based on the behavior of the reactivity variation alone as follows. Initially, when the core is far subcritical, the core power falls as time increases even though the system reactivity is increasing. This happens due to the fact that the system is so far subcritical that it is not able to maintain a steady state neutron population. Then the total system reactivity reaches a point (still subcritical) where the rate of fuel mass increase within the core forces the power to go up (note
85
that the core power is given by an expression which involves a product of the flux and the fuel atom densities within the core. So, even though the flux is still falling (since the core is still subcritical), the rate of increase of the fuel atom densities compensates for this effect and the net result is a power that increases with time). Then at 50 seconds after the initiation of the transition, the system reactivity becomes exactly zero (critical stationary system), and, at this point, within DYNAM, the positive reactivity insertion rate is reduced to a much lower value (see Fig. 3.11). The system power now drops slightly due to fact that, for this circulating fuel system, there is some loss of delayed neutrons while going through the external loop and, even though the stationary system reactivity is zero at time t = 50 seconds, the circulating system reactivity is still negative (system multiplication factor is less than one). Eventually, at t =85 seconds, the circulating system becomes critical and further reactivity addition leads to a rapid power increase.
To understand the behavior of other system parameters (like core temperature, core pressure, etc.) the effect of the balance ofplant also has to be taken into account. The effect of reactivity variation alone cannot explain the behavior of these system variables. For example, Fig. 3.14 shows that the core temperature increases from the very beginning of the transition even though the core power is dropping during the initial phases of the transition. The
1 E+8
_1E+7 S1E+6 S1E+5
1E+4
S1E+3,
IE+2. Q I FlIi
20 40 60 80 100
Time After Startup Initiation (sec)
Figure 312. BPGCR Thermal Power Vs. ime  Case I
Time After Startup Initiation (sec)
Figure 313.
BPGCR Total Fuel (UF6+He) Mass Variation With Time  Case I
I
1400
V
a
~ 1200
U '4
0
Q
U
w
'4 U 1000
0 20 40 60 80 100 120 Time After Startup Initiation (sec)
Figure 314. BPGCR Average Core Temperature Variation with Time  Case I
30.0
20.0
P* 6.OE4 sec
TI= 1.0 sec
0
10.0For the Reactivity Insertion Shown in Fi.31l
o 10.0(Fg
0 20 40 60 80 100 120 Time After Startup Intiation (sec) Figure 315. BPGCR Average Core Pressure Variation with Time  Case I
U 20 40 60 80 100 120
Time After Startup Initiation (sec)
Figure 316. BPGCR Fuel Mass Outflow Rate Variation with Time  Case I
reason for this behavior of core temperature can be directly traced back to the compressor work done on the fuel/working fluid. As mentioned earlier, since this is a Brayton cycle, the compressor work needed to compress the gaseous fuel is considerable. This effect dominates during the initial phase of the transition and the core power (fission power) has very little effect on the core temperature. The work done by the compressor depends on the fuel mass flow rate and, hence, initially, during the rapid insertion of reactivity (addition of fuel to the core), the work done by the compressor is also increasing. Only towards the end of the transition, when the core power becomes significant (above a few megawatts), does the core power starts influencing the behavior of core parameters like core pressure, core temperature, etc.. From Fig. 3.16, we can see that most of the fuel mass flow rate increase occurs during the first 50 seconds (except for the sharp increase towards the end) of the transition. This behavior of the fuel mass flow rate gets reflected in the fuel temperature profilewhere the core temperature pretty much comes to a steady state after the first 50 seconds or so. The discontinuity in the core parameters around 50 seconds after the initiation of the transition can be attributed to the discontinuity in the reactivity insertion rate at this point in time (i.e., a change in the reactivity insertion rate results in a change in the fuel mass inflow rate and, hence, a change in the compressor work output). The core temperature profile does not show any sharp discontinuity at 50
seconds because of the damping that occurs in the external loop (due to the heat exchanger/ radiator).
The type of reactivity insertion described above (pure ramp
insertion) is not a good one to employ during a real transition. As can be seen from Figures 3.12 through 3.16, most of the core power increase (and, hence, the increases in the core temperature, core pressure, etc.) occurs very rapidly during the last few seconds of the transition. This means that some type of drastic controlling action needs to be taken at the end of the transition period to prevent power overshoot. Also, most of the time, the safety features engineered into a reactor are designed to prevent such rapid power increases, and will cause the reactor to be scrammed. So, instead of using just a pure ramp insertion of reactivity, a specific reactivity variation, as shown in Fig. 3.17, was next selected for the purpose of demonstration. This reactivity variation is designed to produce a fast power increase to a predetermined constant power level. As shown in Fig. 3.17, the reactivity variation is a combination of two linear ramps (with different slopes) and an exponentially decaying reactivity variation. By adjusting the peak value of the core reactivity (pna,) and the time at which the exponentially decaying reactivity variation begins (ti), the steady state power level attained can be varied (see Ref.28 for details of this type of reactivity variation. Ref.28 also has details about other types of

PAGE 1
STARTUP AND STABILITY OF A GASEOUS CORE NUCLEAR REACTOR SYSTEM Ey KIRATADAS KUTIKKAD A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY UNIVERSITY OF FLORIDA 1991
PAGE 2
^(9
PAGE 3
ACKNOWLEDGMENTS The author wishes to express his sincere gratitude and appreciation to Dr. Edward Dugan for his patience, constant guidance and encouragement during the course of this investigation. The author has profited immensely by associating professionally with Dr. Dugan. Dr. Dugan's deep understanding of the physics of nuclear reactors and, in particular, his knowledge of reactor kinetics has been a constant source of help to the author during the course of this work. The author also wishes to extend his sincere thanks to the members of his supervisory committee for their willingness to assist in the preparation of this dissertation. Dr. Alan Jacobs needs special mention in this regard. His comments and criticisms were of great help during the final phases of this work. Financial assistance was provided partly by the Nuclear Engineering Department and partly by the Innovative Nuclear Space Power Institute at the University of Florida. The author is thankful to Dr. Nils Diaz, Dr. Ed Dugan, and Prof. Jim Tulenko for this financial support. The author wishes to express his appreciation to a number of his fellow graduate students at the University of Florida. Special thanks go iii
PAGE 4
to Dr. Mathew Panicker and Jerry Welch for the many technical discussions the author had with these two individuals, and to Bob Williams for his cheerful company. Appreciation and affection are extended to the author's mother, sisters, sistersinlaw, brothers and brothersinlaw for their constant support and encouragement. Deepest appreciation also goes to his friendsÂ— whose company made the long time it took to finish this work much more enjoyable, quite memorable and precious. Finally, special thanks are extended to Wale Oladiran, Paul Miller and Charlie McKibben of the University of Missouri Research Reactor for their help and understanding during the final preparation of this manuscript. iv
PAGE 5
TABLE OF CONTENTS Page ACKNOWLEDGEMENTS iii ABSTRACT vii CHAPTERS 1 INTRODUCTION 1 Advantages of Gaseous Fuel Over Solid Fuel 2 Special Advantages of Gas Core Reactors for Space Power Applications 3 Present Work Description 5 Dissertation Organization 9 2 PREVIOUS WORK ON GASEOUS CORE REACTORS 1 1 Previous Gas Core Reactor Work at the University of Florida 14 Present System Description 16 Summary of Previous Work Done in Related Areas 20 3 ANALYSIS OF TRANSITION TO HIGH POWER MODE 23 System Description 25 Selection of External Moderator Thickness 27 Spherical Versus Cylindrical Geometry 31 Estimation of Gas Core Reactor Reactivity Coefficients Using the Equivalent Spherical Geometry 36 Models Employed in the Computer Program 45 Programming Aspects 68 Testing of the Computer Program 70 Results 79 V
PAGE 6
Page 4 BPGCR STABILITY ANALYSIS 104 BPGCR Linear Stability Analysis 135 Inherent Stability of the BPGCR Attached to an External Loop 143 5 CONCLUSIONS, RECOMMENDATIONS AND AREAS OF FURTHER RESEARCH 171 Results 172 Recommendations for Future Work 176 APPENDICES A THE COMPUTER PROGRAM DYNAM 182 B BPGCR SYSTEM EQUATIONS AND THE DERIVATION OF THE BPGCR TRANSFER FUNCTIONS 209 LIST OF REFERENCES 230 BIOGRAPHICAL SKETCH 233 vi
PAGE 7
Abstract of Dissertation Presented to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for the Degree of Doctor of Philosophy STARTUP AND STABILITY OF A GASEOUS CORE NUCLEAR REACTOR SYSTEM By Kiratadas Kutikkad December, 1991 Chairman: Dr. Edward T. Dugan Major Department: Nuclear Engineering Sciences Dynamic aspects of a conceptual Gas Core Reactor (GCR) are investigated in this dissertation. The reactor is a bimodal, closed cycle, circulating fuel system capable of operating in low and high power modes. It uses uranium hexafluoride (UFe) in gaseous form as the fuel, and the fuel gas mixture of UFe and helium also serves as the coolant /working fluid. The gaseous nature of the fuel, plus the fact that this is a circulating fuel reactor, makes the system unique, and the dynamic behavior tends to be quite different from that of conventional, solid core systems. The system behavior while undergoing rapid transition from a low power mode of operation to a high power mode of vii
PAGE 8
operation is investigated with the help of a computer program developed for that purpose. The program, DYNAM, written in standard FORTRAN 77, simultaneously solves the relevant point reactor kinetics, thermod)mamic, heat transfer, and 1D isentropic flow equations of the OCR. A nonlinear stability analysis of the system is also performed using DYNAM. The results of this analysis indicate that the investigated OCR is an inherently stable systemÂ— mainly due to the strong negative fuel mass feedback effect. Due to this strong fuel mass feedback, the response of the OCR to external reactivity perturbations is found to be quite different from that of conventional, solid core reactors. For example, at full power operation, it is observed that the inherent reactor response is a power decrease for a positive, as well as a negative, external step reactivity insertion. To increase GCR power, it is necessary to provide additional positive reactivity during the transient by varying a system operating parameter (for example, by increasing the core inlet pressure). A linear stability analysis of the GCR shows that the linear analysis is often adequate for correctly predicting system stability. However, for certain operating regimes, where the fuel mass feedback coefficient is large, the linear analysis breaks down due to the strong nonlinear effects. For these regimes, the nonlinear analysis has to be performed to reliably predict the system stability. viii
PAGE 9
CHAPTER 1 INTRODUCTION Terrestrial power generation making use of nuclear fission is currently accomplished with solid fueled reactors. The fuel usually is in the form of small pellets made of uranium dioxide (UO 2 ) or in some cases uranium metal or some other compounds of uranium like uranium carbide or nitride. The concept of using nuclear fuel in gaseous form was first suggested in 1955 (Ref.l). Here the fuel can be either uranium plasma (provided the temperature in the system is kept high enough for the plasma to exist) or a compound like uranium hexafluoride (UFe) or uranium tetrafluoride (UF4) in gaseous form. Nuclear fuel in gaseous form has many potential advantages over that in solid form which could eventually make it the fuel form of choice for many special applications including space power. In this dissertation. d)mamic aspects of a conceptual gas core reactor system are investigated. The reactor uses highly enriched UFe or UF 4 in gaseous form as the fuel. The peak fuel temperatures attained are in the range of only a few thousands of degree Kelvin and, hence, the fuel is not in a plasma state. 1
PAGE 10
2 Knowledge of a reactor system's behavior in the nonsteady state is of primary importance in the design and construction of any new type of reactor. The purpose of the research described in this dissertation is twofold. The first goal is to perform the system startup analysis, i.e., analysis and prediction of the system behavior during various startup scenarios and the second goal is to analyze the stability of the system's high power mode of operation, i.e., to establish the conditions under which the reactor can be operated in a safe and reliable manner at high power levels. Advantages of Gaseous Fuel Over Solid Fuel One of the major limitations of conventional solid core pressurized water reactors (PWR) or boiling water reactors (BWR) is their relatively low working fluid temperature. Since the zircaloy cladding of the fuel loses its strength at around 700 K, the working fluid/coolant temperature has to be kept below this value. This is a low value compared to the 1500 K or so reached at the centerline of an average fuel pellet inside the core, and this results in a low power plant energy conversion efficiency compared to a coal or oil fired power plant. The reason for not replacing zircaloy with some other metal or alloy with better mechanical strength is that, except for this drawback, zircaloy is an excellent cladding material with many desirable properties like low neutron absorption cross section, good resistance to corrosion by water at elevated temperatures, and reasonably good resistance to radiation damage. By going to a gaseous
PAGE 11
3 fuel, this serious limitation of solid core systems can be circumvented. A gaseous core system can operate at extremely high temperatures (and high efficiencies) since the fuel, either by itself or when mixed with another gas, can also be the coolant/ working fluid, and, unlike in a solid core system, there may not need to be any heat transfer process from the fuel to a coolant/working fluid. Practical heat removal and materials considerations generally force one to operate the gas core reactor at temperatures in the range of a few thousands of degree Kelvin. An exception to this is the plasma core reactors where the fuel temperatures are in the tens of thousands of degree Kelvin range. Here special arrangements are made to keep the fuel away from the containment wall. A gaseous core system also opens up the possibility of continuous operation of the power plant with online refuelling and online fission product removal. Other advantages associated with gas core reactors include a much less complicated core structure, high fuel utilization, simple fuel management, and inherent safety due to the negative reactivity associated with an expanding gaseous fuel. Special Advantages of Gas Core Reactors for S pace Power Applications The high working fluid temperature plus the gaseous nature of the fuel makes the system quite attractive for space power applications. First of all, a high working fluid temperature implies high efficiency, which in turn implies a compact system for a given power
PAGE 12
4 output. A second advantage is the reduction in the size of the radiator. For any space power source, nuclear or nonnuclear, the problem of efficiently rejecting waste heat is always a major concern. Radiation is the only way of rejecting heat in space, and, for large power systems, the mass of the radiator usually tends to dominate the total system weight (Ref. 2). Thus, there is a strong incentive for reducing the radiator mass. Since the rate of radiative heat transfer varies as the fourth power of temperature, the very high working fluid temperatures possible with gas core reactors offer significant advantages with regard to radiator size and, hence, system mass. A third advantage is the elimination of the accidental criticality concerns. Before the launching of any reactor into space, extreme precautions have to be taken to make sure that the system remains subcritical throughout normal launching operations as well as during any foreseeable accidents during launch. By using a gaseous core system, the fuel in its storage container(s) is (are) in a highly subcritical state during launch. Achievement of criticality requires that the gas be both pressurized and contained or surrounded by an efficient moderator/reflector. Also, since the process of putting a reactor in orbit is an extremely costly affair, one would like the reactor to remain in operation for a long period of time. The gas core reactor offers the possibility of long term operation without building in a large amount of excess reactivity at the time of launch due to its capability for online refuelling and fission product removal.
PAGE 13
Present Work Description The concept analyzed in this dissertation is a gaseous core 5 reactor power system that is designed to operate at two widely different power levels: a low power mode of operation and a burst or high power mode of operation. It is a highly enriched UFe or UF 4 fueled reactor capable of delivering power anywhere from a few MWe up to a few hundreds of MWe. The advantage of using UFe is that it is the only stable compound of uranium that can exist in a gaseous state at a relatively low temperature (sublimation point = 329.5 K at atmospheric pressure). This low melting point implies that the system temperature, including that of the external loop, can be kept at a reasonably low value, and still avoid condensation of the fuel on the inner surfaces of the pipes and ducts. The limitations of UFe are that it is highly corrosive and it undergoes significant dissociation by 2000 K. As for UF 4 , significant dissociation starts only around 5000 K and UF 4 is also much less corrosive than UFe . The problem in using UF 4 as the fuel is that it is very difficult to start the system up because of the high melting and boiling points (1300 K and 1700 K, respectively, at atmospheric pressure) of UF 4 . One possibility is that the system could be started up using UFe as the fuel and then once the fuel temperature gets high enough (above 2000 K), UF 4 naturally becomes the fuel since one of
PAGE 14
6 the major dissociation products of UFe is UF4. A more detailed description of the system is given in Chapter 2. The incentive for designing such a bimodal system primarily came from Strategic Defense Initiative (SDl) requirements. The power requirements specified by the SDIO can be divided into three distinct ranges: a) a station keeping mode requiring a few MWe for a long period of time (710 years). b) an alert or enhanced surveillance mode requiring a few tens of MWe for shorter periods of time (a few hours to a few days) and finally, c) a burst power or defense mode needing a few hundreds of MWe or more for short periods of time (a few tens of minutes up to an hour or so). The bimodal gas core reactor system that is described in this dissertation is expected to meet all of the above power requirements, i.e., the low power mode is capable of meeting the power requirements of both the station keeping phase as well as the alert phase and the high power mode of operation can handle the power requirements during the burst power phase. As a part of the system design, both steady state (static) analysis and dynamic analysis have to be performed. This dissertation is mainly concerned with the dynamic analysis of the system.
PAGE 15
7 Logically, the dynamic analysis can be divided into three distinct phases: a) analysis of the system's transient behavior during the low power mode of operation. b) analysis of the transition from the low to the high power mode of operation, and finally. c) analysis of the system's transient behavior during the high power mode of operation. Dynamic analysis of the low power mode has already been performed at the University of Florida (Ref.3). This work, therefore, concentrates on the remaining two areas. Specifically, the first part of the dissertation looks at the startup of the high power (or the transition to the high power) mode of operation. One of the requirements which any spacebased power system for SDl burst power applications has to satisfy is the capability for rapid startup, i.e.. within a short period of time the system should be able to go safely and in a predictable manner from a low power, stationkeeping mode to a high power mode. SDIO requirements limit this transition time to about 100 seconds, which is extremely short compared to the few hours taken by conventional, terrestrial nuclear plants to reach full power operation starting from low power or hot standby conditions.
PAGE 16
8 The response of the gas core reactor that is analyzed in this dissertation is expected to be quite different from that of a conventional, solid core system. The highly compressible nature of the fuel, along with the fact that this is a circulating fuel system, implies that the fuel mass within the core is not going to be a constant during transients. The interaction between the core pressure, the core temperature, and the fuel mass within the core adds additional nonlinearity to the system equations. It is very important to know the behavior of key system parameters during the rapid transition from the low power mode to the high power mode of operation in order to ensure a safe transition. To analyze the dynamic behavior of the gaseous core nuclear reactor during rapid transitions, a computer program is developed. The program (DYNAM) solves the complete system equations (nonlinear) numerically and predicts the behavior of important system parameters during the transition. It can be a valuable tool for stud)nng the effects of variations in different system parameters on the system response and also for defining safe transition paths. The second part of the dissertation is concerned with the stability analysis of the high power mode of operation, i.e., assuming that a safe transition has already been made from a low power mode of operation to a high power mode of operation, this section looks at the power behavior following any perturbations (primarily, reactivity perturbations) imposed on the steady state, high power system. If the
PAGE 17
9 system is found to be inherently unstable, the system design has to be modified to make it stable. One way to perform the stability analysis is to linearize the nonlinear system equations, and then use the powerful methods of linear stability analysis to determine the system stability. One of the objectives of this section is to determine whether a linear stability analysis is adequate or not for predicting the stability of the gaseous core reactor system. To this end, the results from the linear analysis can be compared against those obtained using DYNAM. If the linear analysis is found to be adequate, one can use it to determine quickly the range of values for reactor physical parameters and equilibrium power levels for which the system is stable, and to predict the bounds on allowable departures from the equilibrium states without having to resort to nonlinear analysisÂ—which can get very involved and time consuming for a system like this. Dissertation Organization Chapter 2 contains a review of the previous work done on gaseous core reactor power systems. It includes a review of the gas core reactor work that has been going on at the University of Florida for the past 15 years or so. It also includes a detailed description of the bimodal reactor that is the topic of this dissertation. This is followed by a chapter which describes the transition or startup analysis of the high power system. The models used in the computer code developed for analyzing the transition, as well as some
PAGE 18
10 of the limitations of the code, are discussed. This chapter also presents some of the representative results generated by the code. Chapter 4 contains the theory and a description of the methods used for the nonlinear and linear stability analysis of the burst power system. It also contains the results of the stability analysis. The last chapter presents the conclusions. It also contains some suggestions for future work, and modifications that can be made to the computer program, DYNAM, written for analyzing the gas core reactor.
PAGE 19
CHAPTER 2 PREVIOUS WORK ON GASEOUS CORE REACTORS Research on gaseous core reactors has included many independent, conceptual phase efforts for the past 30 years or so. During this time, a wide variety of system configurations and physical arrangements for the reactor core, energy conversion schemes, and fuel/working fluid form have been suggested, and a wide range of scientific feasibility studies, mostly analytical and a few experimental, have been conducted with encouraging results. As mentioned in Chapter 1, in the United States the gas core reactor concept was first suggested by Bell in 1955 (Ref. 1). He performed preliminary static neutronics calculations on a spherical core fueled with UFe gas and surrounded by a moderator/reflector made of either heavy water (D 2 O), beryllium or graphite. The impetus for much of the early work on gas core reactors came from the National Aeronautics and Space Administration (NASA). Primarily, NASA was interested in the possibility of applying such systems for rocket propulsion in space. A number of concepts were evaluated by NASA and out of these, only two were selected for further study (Ref.2). These two were the opencycle coaxial flow system and the closedcycle nuclear light bulb concept. The general idea in both these concepts is to use the energy
PAGE 20
12 from a fissioning plasma (of partially ionized uranium with average fuel temperatures in excess of 10,000 K) to radiatively heat the propellant which is then exhausted through a nozzle to generate the required thrust. The original work on the opencycle coaxial flow reactor was performed by Rom (Ref. 4) and Ragsdale (Ref. 5) under the direction of the NASA Lewis Research Center. The coaxial flow concept involves a spherical cavity with a low velocity inner stream of fissioning uranium metal plasma. This hot plasma is kept away from the cavity walls by a high velocity propellant gas (typically hydrogen) stream. The hot inner fuel plasma transfers thermal energy to the outer propellant stream by both convection and radiation heat transport. In this design, the heat transfer occurs essentially unimpeded, but the drawback is that the propellant stream can physically mix with the fuel stream leading to a loss of nuclear fuel, A more detailed description of the system can be found in Schneider and Thom (Ref. 6). To avoid the mixing of the propellant and fuel streams, in the nuclear light bulb concept, the fissioning plasma is confined within a transparent (generally quartz) wall. The fuel is kept away from this transparent wall by injecting a buffer gas (argon, neon, etc.) tangentially into the cavity. Here, the energy transfer from the plasma to the propellant is by radiation heat transfer alone. Major work on this concept was performed at the United Aircraft Research Laboratory (UARL). Lantham (Ref, 7) performed a series of calculations
PAGE 21
13 for the nuclear light bulb reactor rocket engine under the auspices of UARL. The application of gas core reactors to terrestrial electric power generation also had some early investigations. In 1969, Gritton and Pinkel (Ref. 8) of the Rand corporation reported one such study in which they evaluated a 4000 MWt spherical gaseous core system. The central gas core chamber was surrounded by a moderating/reflecting region and by banks of energy conversion devices. In another study, Williams and Shelton of the Georgia Institute of Technology (Ref,9) looked into the possibility of coupling gas core reactors with a Magneto Hydrod 3 mamic (MHD) power conversion system for terrestrial power generation. They studied the possibility of coupling both the coaxial flow system and the nuclear light bulb system to an MHD generator. They concluded that gas core reactors, in conjunction with MHD generators, are capable of generating electricity in a safe and efficient manner. Up until then, the main obstacle preventing MHD generators from achieving their full potential was the low electrical conductivity of the working fluid. The study conducted at the Georgia Institute of Technology concluded that gas core reactors are capable of providing the high temperatures and the nonequilibrium ionization necessary for high electrical conductivities and, hence, for high MHD efficiencies.
PAGE 22
14 Previous Gas Core Reactor Work At The University of Florida As mentioned above, the major early gas core reactor research efforts in the United States were closely related to the NASA's space power needs and the program was essentially abandoned in the early 1970s. But, at the University of Florida (U of F), research on gas core reactors continued without much interruption during the 1970s and the 1980s. During these years many different concepts were studied and characterized at the U of F. The work on gaseous core reactors started at the U of F in the early 1970s when Schneider and Ohanian (Ref. 10) proposed a nonsteady reactor called the Pulsed Nuclear Piston (PNP) engine. This was a pulsed system operating on a thermodynamic cycle very similar to the internal combustion engine with the fissionable gas as the working fluid. This piston engine consists of a cylindrical core enclosed by a cylinderpiston arrangement just like in an internal combustion engine. This is surrounded by a moderating/reflecting region. A mixture of UFe and helium is used as the fuel and the working fluid (the helium is added to enhance the thermodynamic and heat transfer characteristics of the primary working fluid). The fuel is introduced into the core, and is then compressed by the piston to form a supercritical system. Fission energy is rapidly released within the cylinder thereby increasing the temperature and pressure of the fuel
PAGE 23
15 considerably. The fuel is then exhausted from the cylinder and the cycle is repeated. Neutronically, the core goes from a highly subcritical state to a critical state and then to a highly supercritical state. The energy released by the fissioning gas can be extracted both as mechanical power (from the motion of the piston) and as heat from the circulating gas (e.g., through a heat exchangerturbine combination). Because of this double mode of energy extraction, the system efficiency can be as high as 50% (Refill). An extensive neutronics and energetics analysis of this PNP concept was performed by Dugan (Ref. 12). A pulsed system which is only slightly different from the PNP was also studied by Diaz et al. (Refill) in 1979. This is the Pulsed Gas Generator (PGG) concept. The PGG is very similar to the nuclear piston concept except for the elimination of the piston and making the core volume fixed. The fuel (which is again a mixture of UFe and helium) is cyclically introduced into the core just like in a PNP. The pulsing effect is now achieved with the help of rotating absorber rods placed in the external moderating region. The PGG concept has the advantage of mechanical simplicity since it requires no moving piston. However, the efficiency attained by the system is lower than the PNP due to the fact that only heat energy (turbine power) is extracted from the system now. A more detailed description of the system can be found in Diaz et al. (Refill).
PAGE 24
16 A gaseous core system significantly different from the pulsed systems (PNP and PGG) was next in line to be analyzed at the U of F. This is the Heterogeneous Gas Core Reactor (HGCR). This novel concept in gas core reactor design was suggested in 1974 by Diaz and Dugan (Ref. 13). A heterogeneous core is formed by arranging bundles of moderator/ coolant channels in a fissionable gaseous core region. Thus the HGCR has 'internal moderation' as opposed to the external moderation in the PNP and PGG. The core heterogeneity leads to significant thermal hydraulic and energetic advantages along with improved technological feasibility compared to other gas core reactor concepts. Ki Han (Ref. 14) performed static neutronic analysis and some kinetic analysis of this system. Present System Description One of the latest in the series of gas core reactor concepts that has been analyzed at the U of F is a bimodal nuclear power system. A schematic of the core proper is shown in Fig. 2.1. Typical reactor dimensions and other key system parameters are given in Table 2. 1 . The development of the system design and characterization are being done under the auspices of the Innovative Nuclear Space Power Institute (INSPI) at the U of F. As shown in Fig. 2.1, the system consists of two identical halves. Basically, each half of the system consists of a central fissioning region (called a Burst Power Gas Core Reactor or BPGCR) surrounded by a ring of mini cylinders (called Pulsed Gas Core Reactors or PGCRs)
PAGE 25
17 Figure 21. Bimodal Gas Core Reactor Schematic
PAGE 26
Table 21 . Bimodal Reactor Typical Dimensions and Operating Conditions. PGCR Data (Each Half): Number of PGCFte 8 12 Height (m) 1.0 Diameter (m) 0.25 Thermal Power/PGCR (MWt) BPGCR Data (Each Half): 3.0 Diameter (m) 1.0 Height (m) 2.0 Thermal Power (MWt) 200.0 External Beryllium Thickness (m) 0.5 Fuel (UFe + He) Mass (Kg) 10.0 Average Core Pressure (atm) 20 30 Average Core Temperature (K) 1500.0 Helium Mole Fraction in the Fuel 0.93 U 235 Enrichment in the Fuel 0.85
PAGE 27
19 with a region of moderating beryllium or beryllium oxide in between. The PGCRs are also surrounded by a reflecting and moderating layer of beryllium or beryllium oxide. Cooling channels, carrying some inert gas like helium, are provided in the moderator region to keep its temperature below the melting point. Each half of the core exhausts fuel gas into a central duct and nozzle configuration which then directs the fuel gas into a disc MHD generator. The core fuel gas temperature is in the range of 20004000 K in this concept. This high temperature of the fuel, along with the significant amount of fission fragment induced ionization taking place in the MHD duct, is expected to provide adequate electrical conductivity and reasonable MHD generator efficiencies. The PGCRs are identical to the pulsed gas generators described earlier. Each PGCR consists of a cylindrical core of roughly a meter in height and 20 to 30 centimeters in diameter. The core is embedded in a beryllium or beryllium oxide moderating/reflecting region. Each half of the bimodal power system includes from 8 to 12 PGCRs, each of which is capable of delivering roughly 3 MWs of thermal power (Ref. 15). The low power mode of operation of the bimodal system is defined, for the purpose of this dissertation, as the case where only the PGCRs are suppl 3 dng power. The central chamber is left empty or filled with an inert gas. The power from the PGCRs is expected to be sufficient to meet the power requirements of both the station keeping mode and the alert mode as mentioned in Chapter 1. Panicker (Ref.3)
PAGE 28
20 has performed extensive static and dynamic neutronic analysis of the PGCR operation for this bimodal power system. Using Monte Carlo techniques, he has also studied the coupling effects among the PGCRs and between the PGCRs and the BPGCR. The BPGCR is a cylindrical chamber about 2 meters in height and roughly a meter in diameter (each half). This core is surrounded by a beryllium or beryllium oxide region which acts as an external moderator and reflector. Each half of the BPGCR is capable of delivering 100s of MWt of power, and during the high power (or burst power) mode of operation of the bimodal system, the BPGCR will be supplying most of the system power. The focus of this dissertation is the dynamic analysis of the BPGCR. It should be noted that since this system is only in a conceptual design state, variations of the presented design are also being examined. However, for the purpose of this dissertation, this base design is selected. The analysis that is performed in this dissertation should also be valid for other designs, provided the concept is not changed drastically. Summary of Previous Work Done in Related Areas As mentioned earlier, a systematic study of gas core reactors was never undertaken in the United States, As a result of this, plus the fact that much of the work related to gas core reactors (spacebased systems in particular) was discontinued in the early 1970s, only a couple of studies on the startup and stability of gas core space nuclear
PAGE 29
21 power systems are available. Although there are many such studies for conventional, terrestrial nuclear plants, they cannot be extended directly to the gas core space nuclear systems because of the different designs and quite different demands on these systems. Gresho et al. (Ref. 16) and Birken (Ref. 17) deal with the startup analysis of the SNAP reactor systems (SNAP 2 and SNAP 8). These were small, solid core thermal reactors developed for meeting the power requirements of various space missions. They were capable of generating only a few kilowatts of electrical power. Since these systems were not meant for any defense applications, their startup requirements were much less severe than those of the bimodal gas core reactor system that is analyzed in this dissertation. For example, the SNAP 2 system took about 6 hrs. in going from a far subcritical state to a full power, steady state operation. This long startup time allowed the transition to be made in a much more controllable manner than for the bimodal gas core system. For the startup analysis of these systems, mostly analog computer techniques were used, and the only feedback effects considered were the fuel and grid plate temperature coefficients. The startup analysis of the gas core system is expected to be very different from the SNAP system analysis because of factors like the gaseous nature of the fuel, the excore circulation of the fuel gas mixture (which is also the coolant/ working fluid), delayed neutron losses due to the fuel circulation, etc.
PAGE 30
22 The SNAP10 reactor, which belongs to the same class of reactors as the SNAP2, was launched into space in 1965. The reactor was started up while in orbit, and it operated successfully in the space environment for more than a month. Actual experimental data are available on the startup and shutdown of this reactor. A more detailed description of the system can be found in Angelo and Buden (Ref. 2). In the startup analysis of the nuclear light bulb engine, the authors (Ref. 18) did consider rapid startup situations (60, 300 and 600 sec. power ramps). Since the objectives of the study were to estimate the thermal stresses involved, the form of fuel best suited for the system, etc., the neutronic analysis of the startup was performed using only simplistic models. The analysis neglected delayed neutrons completely and no fuel circulation effects were taken into consideration. Also, they looked at only linear power ramps. Lantham et al. (Ref. 19) studied the dynamic response of the nuclear light bulb engine to perturbations induced while the system is operating at full power. To obtain the system response, they solved the full set of coupled thermal, fluid dynamics and neutron kinetics equations. Ki Han (Ref. 14) performed some stability analysis of the high power heterogeneous gas core reactor. His analysis was preliminary in nature and he treated fuel circulation effects in a very crude way. He also did not include the balanceofplant in his analysis.
PAGE 31
CHAPTER 3 ANALYSIS OF TRANSITION TO HIGH POWER MODE Introduction The purpose of the work described in this chapter is to analyze the transition of the bimodal power system from a low power, stationkeeping mode of operation to a high power, burst mode. As mentioned in Chapter 1, this transition has to be made very fastÂ— in about 100 secondsÂ— and it is very important that key system parameters like core power, core temperature, core pressure, fuel mass flow rates, etc. be closely controlled during the transition to make sure that no maximum limits on any of these variables are exceeded. The information obtained as a result of this analysis can also be used in the overall system design. For example, the external beryllium (or beryllium oxide) moderator/reflector has to be continuously cooled to remove heat deposited in it due to convective, radiative, and neutron/gamma heating from the core so as to maintain this region at the desired temperature. The heat deposition in the moderator needs to be known to establish the required coolant flow rates through the beryllium to prevent it from overheating and 23
PAGE 32
24 melting. The temperature change of beryllium also introduces feedback reactivity which needs to be taken into account when determining the system's neutron multiplication factor. For the purpose of analyzing the transition, a computer program was developed. This program can accept any userspecified reactivity variations as input, and predict the behavior of the core parameters of interest during the transition. The program is written in standard FORTRAN 77 and can run either on PCs or mainframes. In practice, the userspecified input reactivity variation during the transition can be accomplished in a number of ways. For example, as in a conventional terrestrial nuclear power plant, the fuel mass inside the core can be kept fixed by maintaining constant and equal fuel mass inflow and outflow rates and then with the help of rotating absorber drums in the external moderator/reflector, the necessary reactivity variation can be brought about. Another way of varying the core reactivity is to vary the fuel mass inside the core. This is possible because of the gaseous nature of the fuel, plus the fact that this is a circulating fuel reactor. The fuel mass variation inside the core can be achieved in a number of different ways. One wayÂ— which is described later in this chapter to demonstrate the working of the computer programÂ— is to control the fuel mass inflow rate into the core by controlling the fuel inlet pressure to the core. The fuel mass outflow rate from the core is then determined by the conditions existing within the core (temperature, pressure, etc.) and by the downstream
PAGE 33
25 conditions. The code can be modified, with minimum effort, to handle any other plausible startup scenarios. System Description The schematic diagram of the anal)^ed reactor core is given in Fig. 2.1. Typical dimensions of the core are given in Table 2.1. For the purpose of this dissertation, it is assumed that the PGCR's are voided (unfuelled) during the transition to high power. A schematic diagram of the complete system, including the balanceofplant, is shown in Fig. 3.1. In Fig. 3.1, the two halves of the BPGCR core are represented as a single core. The output from the core enters the MHD duct after passing through a convergingdiverging nozzle (not shown in Fig. 3.1). The purpose of this nozzle is to accelerate the fuel exiting from the core to supersonic velocities before entering the MHD generator (high working fluid Mach numbers are necessary for reasonable MHD generator efficiencies). From the MHD duct, the fuel passes through a supersonic diffuser (where its velocity gets reduced considerably) and a radiator heat exchanger before getting pumped back into the core. A fuel reservoir is located before the compressor to provide any excess fuel needed while starting up the system, and to act like a sink to retain any excess fuel from the loop during burst power operational transients or shutdown. The heat exchanger transfers heat to a radiator for rejecting it to space.
PAGE 34
26 Figure 31. Schematic of Complete BPGCR Power System Including BalanceofPlant.
PAGE 35
27 Referring to Fig, 2.1, it can be seen that the reactor consists of two identical halves, and, since the PGCRs are assumed to be shut down during the transition to high power, the reactor can be essentially treated as two large identical cylindrical BPGCR's. The analysis is simplified by analyzing the behavior of only onehalf of the system during the transition. The Justification for this simplification is the assumption that the neutronic coupling between the two cores is expected to be small (this coupling actually depends largely on the separation distance between the two halves, i.e., on the thickness of the separating beryllium or beryllium oxide region). The core power and mass flow rates obtained from this analysis can then be multiplied by a factor of two to get the corresponding values for the entire system. Selection of External Moderator Thickness A 'cavity reactor' configuration is one where the gaseous fuel region is surrounded by a moderating/ reflecting region and, unlike in a conventional solid core reactor, almost all of the neutron moderation takes place in this region external to the fuel region. The BPGCR is such a cavity reactor with a beryllium (or beryllium oxide) moderator surrounding the core. An adequate thickness of the moderating layer is crucial for achieving criticality. If this thickness is too small, the fission neutrons created within the core do not have a good chance of being moderated in the beryllium region and reflected back into the core to cause further fissions. Most of them will leak out of the system.
PAGE 36
28 and, if the thickness is smcill enough, adequate neutron thermalization also will not occur. So, to maintain criticality, either the fuel mass inside the core has to be increased, or the core size has to be increased. At the other extreme, the thickness of the beryllium region cannot be made arbitrarily high due to system weight considerations. This is especially true in space power applications where one of the main goals of the system design is to keep the system weight down. So it is important to find an optimum thickness for the beryllium moderator. To find the optimum beryllium thickness, the cylindrical core is first mocked up in spherical geometry by conserving the total core volume (see the next section for a discussion on the validity of this procedure). Then four group, onedimensional Sn transport theory keff calculations using XSDRNPM (Ref. 20) are performed for various beryllium thicknesses. The obtained results are given in Fig. 3.2. These calculations are performed for an assumed core pressure of 50 atmospheres and a core temperature of 3000 K. The fuel atom densities can be estimated by assuming ideal gas behavior for the UFeHe mixture. The calculated densities are given in Table 3.1. As can be seen from Fig. 3.2, after the beryllium thickness exceeds about 50 cm., no significant gain in system keff is obtained by increasing the moderator thickness further. So, for the purpose of this analysis, a 50 cm. beryllium thickness is selected.
PAGE 37
Table 31. Core Fuel (UFeHe) Atom Densities Nuclide Atom Density (#/cc) U235 7.276 E+18 U238 1.283 E+18 F 5.135 E+19 He 1.137 E+20 Helium Mole Fraction in the Fuel = 0.93 Uranium Enrichment = 85 atom% Average Core Pressure = 50 atm Average Core Temperature = 3000 K
PAGE 38
1.20 30 JJ33I HOOda Figure 32. BPGCR keff Vs. Beryllium Moderator/Reflector Thickness
PAGE 39
31 Spherical Versus Cylindrical Geometry Most of the keff calculations in this dissertation were performed by mocking up the actual cylindrical geometry with an 'equivalent' spherical geometry, and then performing onedimensional, multigroup discrete ordinates calculation using XSDRNPM. The 'equivalent' spherical geometry is formed by conserving the total core volume and the thickness of the external beryllium moderator/reflector. The reason for using a spherical geometry instead of the actual cylindrical geometry is that it is very time consuming and, hence, quite expensive to perform a 2D or 3D cylindrical geometry keff calculation on an externallymoderated gas core reactor using any standard neutronics code. One way of getting around this problem is by performing a 1D (radial dependence only) cylindrical or spherical geometry calculation. When such 1D calculations are perfoirned with cylindrical geometries, a transverse buckling correction factor is generally employed to account for the perpendicular (Zdirection) leakage. This buckling correction procedureÂ— which is routinely used in solid core reactor calculationsÂ— is of no value with gaseous core reactors because of the fact that the flux within the core of a gaseous core reactor is generally very flat. A spherical geometry mockup provides a better onedimensional approximation for the cylindrical gas core reactors since the problem of perpendicular leakage is absent for this configuration.
PAGE 40
32 To check the effectiveness of the equivalent spherical geometry mockup in predicating the kgff of the actual cylindrical core, a series of calculations were performed using XSDRNPM and MCNP (Ref.21). The primary tool used for keff calculations in this dissertation is XSDRNPM. It is a one dimensional, discrete ordinates transport theory code from Oak Ridge National Laboratory which can perform multigroup kgff calculations. MCNP is a Monte Carlo code which can perform continuous and discrete energy neutron, photon, and coupled neutron/photon transport through arbitrary 3D geometries. The MCNP calculations are taken as the standard, and the XSDRNPM results are compared against this. The XSDRNPM calculations were performed using the IBM 3090400 mainframe computer at the University of Florida and the MCNP calculations were performed using the CRAY XMP/ 48 Supercomputer at the San Diego Supercomputing Center. First, a set of keff calculations were performed using XSDRNPM and MCNP on a simple, two region spherical system. The system consists of a spherical core of about 72 cm. radius and a beryllium reflector thickness of 50 cm. For this set of calculations, both XSDRNPM and MCNP treat the same spherical geometry. The results of these calculation are given in Table 3.2. An identical geometry comparison was undertaken to eliminate any geometry effects, and to pinpoint sources of error arising from neutron cross section library
PAGE 41
33 Table 32. MCNP keff Versus 123Group, S4P3 XSDRNPM kefr Geometry Calculation Type keff % Difference Spherical Continuous Energy MCNP 1 1.1140 1.01% Spherical 123 Group XSDRNPM 2 1.1027 1. Estimated Relative Error for MCNP Result = 0.7% 2. Flux and keff Convergence Criterion for XSDRNPM Calculation = 1.0 E4 differences, including energy structure differences (see discussion below). Calculations performed at the University of Florida, using XSDRNPM, use a neutron cross section library which comprises data taken primarily from ENDF/B111 and ENDF/BfV compilations whereas MCNP uses mostly ENDF/BV data maintained at the San Diego Supercomputing Center. The MCNP code is capable of performing continuous energy neutron transport whereas XSDRNPM can perform only multigroup neutron transport. With the cross section library maintained at the U of F, the finest group structure available is 123 neutron groups. The 123group energy structure is quite accurate, and is almost as good as
PAGE 42
34 a continuous energy case. So, the results given in Table 3.2 indicate that the error due to cross section library difference is quite small since the 1% or so difference shown in this table is within the statistical fluctuations of the MCNP results. The next set of calculations were performed to get an estimate of the error involved in using a 4group energy structure compared to a finer structure, for e.g., a 26group or a 123group. Table 3.3 shows the difference in keff values obtained from a 123group, a 26group, and a 4group XSDRNPM calculations. The table also gives the CPU time taken for these three cases on an IBM 3090 mainframe. The geometry and the cross section library used are the same for all the three cases. As can be seen from the table, the 4group calculations predicts a keff value which is 3,1% higher than a 26group calculation and the 26group calculation predicts a keff which is only 0.04% higher than a 123group calculation. Thus, it can be seen that the error involved in using a 4group energy structure is small enough, and the computer time savings involved is significant enough to justify the use of a 4group energy structure. Finally, 'equivalent' spherical geometry XSDRNPM keff calculations (both 123group and 4group) were compared against a cylindrical geometry continuous energy MCNP calculation. The 'equivalent' spherical geometry was formed by conserving the total core volume and keeping the external beryllium thicknesses the same.
PAGE 43
35 Table 33. XSDRNPM S 4 P 3 keff Calculations 4Group Versus 26Group Versus 123Group. Geometry Calculation TVP^ keff % Difference 1 CPU Time (sec) Spherical 123 Group 1.1027 311.48 Spherical 26 Group 1.1031 0.04% 21.01 Spherical 4 Group 1.1370 3.10% 0.53 1. CPU Time on an IBM 3090400 Mainframe for convergence Levels of 1.0 E4. Table 34. Spherical Geometry XSDRNPM S 4 P 3 keff Versus Cylindrical Geometry MCNP keff. Geometry Calculation Type keff % Difference Spherical 123 Group XSDRNPM 1.1027 1 . 2 % Spherical Continuous Energy MCNP 1.0893 Cylindrical 4 Group XSDRNPM 1.1370 4.4% Estimated Relative Error for MCNP Result = 0.7% Flux and keff Convergence Criterion for XSDRNPM Calculation =1.0 E4.
PAGE 44
36 The cylindrical core is 2 meters in height and a meter in diameter with a 50 cm. of beryllium surrounding it. The 'equivalent' spherical system has a radius of 72. 1 cm. and an external beryllium thickness of 50 cm. The results of the keff calculations are given in Table 3.4. There is a 4.4% difference between the MCNP result and a 4group XSDRNPM result. From the previous two comparisons (Tables 3.2 and 3.3), we know that for the same geometry, a 4group calculation predicts a keff value which is about 3% higher than a 123group calculation. In addition, we also know that the cross section library differences are responsible for small differences between the keff values predicted by XSDRNPM and those predicted by MCNP for the same geometry. So, it can be concluded that the error due to geometry difference (cylinder vs sphere) in this particular case is around 2%, This error should be smallest for "square" cylinders (i.e., cylinders with their height equal to their diameter), and should increase as cylinders become elongated. Estimation of Gas Core Reactor Reactivity Coefficients Using the Equivalent Spherical Geometry Reactivity coefficients are very important neutronics parameters of a nuclear reactor, and, when designing any new type of reactor, these coefficientsÂ— which are specific to that systemÂ— need to be determined carefully. They are used in the stability analysis of the system, and the reactor can be stable or unstable at certain power levels depending on the magnitude and sign of these coefficients.
PAGE 45
37 Values for these coefficients are also needed for the startup analysis to determine the total system reactivity during startup. Three major reactivity coefficients associated with the BPGCR are the fuel density coefficient (or equivalently, the pressure coefficient) of reactivity, the moderator temperature, and the fuel temperature coefficients of reactivity. Using the selected equivalent spherical geometry for the gas core reactor, values for these three feedback reactivity coefficients were estimated. The results are given in Tables 3.5 through 3.8. The fuel density coefficient is the most important reactivity coefficient for gas core reactorsÂ— where the fuel is in the form of a compressible fluid. In conventional reactors (FWRs or BWRs), where the fuel is in solid form, this reactivity coefficient is of much less importance because of its much smaller magnitude. The estimated values of the fuel density coefficient of reactivity for the gas core reactor are given in Table 3.5. It is a prompt and positive coefficient. But it should be noted that, as a function of power or temperature, this positive coefficient translates Into a negative effect, i.e., if the core power or fuel temperature increases, the density of the fuel within the core is likely to go down (since the system is an open flow system, an increase in temperature need not necessarily imply a decrease in the core fuel density) and, hence, the system reactivity will also go down. The second reactivity coefficient associated with the fuel inside the core is the fuel (Doppler) temperature coefficient of reactivity.
PAGE 46
38 Table 35, Fuel (UFeHe) Density Coefficient of Reactivity Fuel Density (gm/cc) 1 keff Fractional Change in keff per Unit Density Change (Ak/k per gm/cc) 1,0 0,3752 + 450,50 2,0 0,5934 + 169,50 4,0 0,8356 + 72,79 6,0 0,9668 + 40,78 8,0 1,0490 + 25,50 10,0 1,1039 1. 4Group XSDRNPM S 4 P 3 Calculation on a Spherical Core of Radius 72,1 cm, with a Beryllium Thickness of 50 cm.
PAGE 47
39 Table 3.6 shows the system kefr as a function of the fuel temperature. The fractional change in keff per degree change in the fuel temperature is also given in this table. This fuel temperature coefficient of reactivity is a prompt coefficient (i.e., it follows the changes in the fuel temperature almost instantaneously) and. hence, is a very important coefficient in conventional solid core reactors from a safety standpoint. Table 36. Fuel (UFeHe Mixture) Temperature Coefficient of Reactivity Fuel Temperature (K) 1 keff Fractional Change in keff per Degree K (Ak/k per degree K) 500 1.09867 1.766 E 6 1000 1.09770 +4.191 E7 1500 1.09793 6.558 E7 2000 1.09757 1. 123Group XSDRNPM S 4 P 3 Calculation on a Spherical Core of Radius 72.1 cm. with a Be Moderator Thickness of 50 cm. keff and Flux Convergence Level=1.0 E5. Increasing the fuel temperature affects the system reactivity due to changes in the fuel's microscopic cross section; this is a Doppler
PAGE 48
40 effect in the fuelwhere resonance interaction increases due to the broadening of the resonances. In conventional reactors, the fuel enrichment is very low (34% in U235) compared to the 8085% enrichment used in the gas core reactors described in this dissertation. With such low enrichment, the fuel temperature coefficient of reactivity is dictated mainly by the Doppler effects in the huge capture resonances of the U238 nuclide. So, at all temperatures, the fuel temperature coefficient is negative. In the case of the gas core reactor considered in this dissertation (where the fuel is highly enriched in U235), the fuel temperature coefficient is dictated by very complex phenomena, where the Doppler effect in the capture resonances of U238 and U235 compete with that in the fission resonances of U235. Apart from this, there is a neutron flux spectral effect alsoÂ— where the increased resonance interactions tend to change the neutron flux spectral distribution within the core. Some efforts were made to pinpoint the source of the oscillatory behavior of the fuel temperature coefficient of reactivity (see Table 3.6). But since no definite conclusion could be drawn from the limited study undertaken, this task is left as a suggestion for future work. A change in the moderator temperature affects the system reactivity through two separate mechanisms. The first one is the 'spectral effectÂ’. Since almost all of the neutron thermalization is taking place in the beryllium moderator, an increase in the moderator temperature results in a harder thermal neutron spectrum within the
PAGE 49
41 core (due to the increased vibrational energy of the moderator atoms). A harder thermal neutron spectrum in turn implies a lower overall thermal neutron absorption (fission plus radiative capture) probability. This part of the beryllium moderator temperature coefficient of reactivity is given in Table 3.7. The values given in this table were obtained by repeating keff calculations in which only the beryllium thermal cross section data were changed, i.e., basically using beryllium cross sections with the thermal scattering kernels processed at different temperatures. Changing the thermal scattering kernel of the beryllium will alter the thermal neutron spectrum within the core and, hence, will change the neutron multiplication factor. All other parameters, including the beryllium density, were kept constant for these calculations. As can be seen from Table 3.7, the spectral component of the beryllium temperature coefficient of reactivity switches from a small positive value at low temperatures to an even smaller, but negative, value above 600 K, and then becomes larger and more negative with increasing moderator temperature. From a safety point of view, it is desirable to have this coefficient negative throughout the entire temperature range of interest. The positive coefficient at low temperature is not much of a concern because of the fact that the temperature range where the coefficient is positive is quite small, and at normal operating power levels of the system, the beryllium temperature is expected to be significantly higher than 600 K. Also, in
PAGE 50
42 practice, due to the large thermal inertia associated with the beryllium (or beryllium oxide) moderator/reflector, this reactivity coefficient has a very long time constant associated with it and, hence, will contribute a negligible amount of feedback reactivity during any fast transient. Table 37. Spectral Component of Beryllium Moderator Temperature Coefficient of Reactivity. Beryllium Physical Temperature (K) keff^ Fractional Change in keff per Degree K (Ak/k per degree K) 296 1.08638 +3.485 E5 600 1.09789 2.186 E 6 800 1.09741 1.549 E5 1000 1.09401 1. 123 Group XSDRNPM S 4 P 3 Calculation on a Spherical Core of Radius 72.1 cm. with a Beryllium Moderator Thickness of 50 cm. keff and Flux Convergence Level =1.0 E5. The second mechanism through which the moderator temperature affects the system reactivity is the density effect, i.e., the decrease in the physical density of the moderator with increasing temperature. The beryllium moderator density coefficient of reactivity is given in Table 3.8. This coefficient has a negative effect with respect
PAGE 51
43 to temperature (i.e,, as the physical temperature of the moderator increases, the beryllium density will go down and, hence, the keff will also decrease). Table 38. Beryllium Moderator Density Coefficient of Reactivity. Beryllium Physical Temperature (K) Beryllium Density (gm/cc) keff^ Fractional Change in keff per Degree K (Ak/k per degree K) 296 1.8477 1.08994 7.501 E 6 600 1.8249 1.08743 9.196 E 6 800 1.8063 1.08543 1.230 E5 1000 1.7863 1.08276 1. 123Group XSDRNPM S 4 P 3 Calculation on a Spherical Core of Radius 72. 1 cm. with a Beryllium Moderator Thickness of 50 cm. Flux and keff Convergence Level =1.0 E5. For obtaining the moderator and the fuel temperature coefficients of reactivity, 123group XSDRNPM keff calculations were performed. For the fuel density coefficient of reactivity only 4group calculations were done. The reason for selecting a finer (and, hence, a more accurate) 123group calculation for the moderator and the fuel temperature coefficient determination is that the magnitudes of these
PAGE 52
44 coefficients are quite small (compared to the fuel density coefficient of reactivity), and if a crude energy structure were to be used for calculating them, it may not be able to pick up the energy sensitive effects accurately enough. Also, for the same reason, very stringent flux and keff convergence criteria (105 or less) were imposed for the calculation of these two coefficients. Compared to the moderator and fuel temperature coefficients, the density coefficient of reactivity is quite large (compare Table 3.5 with Tables 3.6, 3.7, and 3.8) and, hence, only 4group calculations were performed to obtain this coefficient. For the transition to burst power analysis, only the fuel density coefficient of reactivity is taken into account. The fuel temperature coefficient is neglected because of its small magnitude compared to the fuel density coefficient, and the beryllium temperature coefficient is neglected during the transition because of the long time it takes for the beryllium to undergo a significant temperature change. The change in beryllium temperature during the transition to burst power was calculated for various scenarios and is given later in this chapter. It was found that, even for the case where no heat is removed from the beryllium moderator, the moderator temperature change during the 100 sec. transition is not significant enough to warrant the inclusion of the beryllium temperature coefficient of reactivity in the transition. In contrast, when the system is operating at high or burst power levels, and full power transients are analyzed, the moderator temperature
PAGE 53
45 coefficient is one of the significantÂ— although slowÂ— feedback effects. So, in Chapter 4, where the stability analysis of the system is described, this coefficient is included. Models Employed in the Computer Program Modelling of the Core Geometry As mentioned earlier, the PGCRs are assumed to be shutdown during the transition to high power. This makes it possible to treat the reactor as two identical halves (BPGCRs). The analysis is further simplified by assuming a weak neutronic coupling between the two cores and, hence, only one half of the system is treated during the transition. One half of the BPGCR is shown in Fig. 3.3. The core is taken to be 2 meters in height and a meter in diameter. It is surrounded by a 50 cm. thick beryllium moderator /reflector (not shown in the figure). As shown in Fig. 3.3, the BPGCR is modelled as an 'open flow' system with no valves at the core inlet or exit; in a real system, there may be some sort of valving or orificing at the core inlet for controlling the fuel flow into the system and, if necessary, for shutting the system down. Similarly, there might be a valve or orifice at the core exit also to control the fuel flow out of the core. The flow characteristics will be altered drastically if the flow through the core is controlled in this manner and, hence, the fluid flow analysis performed here would have to be modified to take this into account. A converging nozzle is provided at the core inlet. At the core exit there is a
PAGE 54
46
PAGE 55
47 convergingdiverging (CD) nozzle. The converging nozzle at the core inlet is expected to provide an 'inherent' control mechanism while the system is operating at full power. For example, if the system power level goes up for some reason, the core pressure and temperature will also go up, and this will cause a reduction in the fuel inflow rate (provided the flow at the inlet is nonchoked and the core inlet pressure is kept constant). The reverse action occurs if the system power drops for some reason. The CD nozzle at the core exit also will provide a similar controlling action. The exit nozzle is assumed to be operating under choked conditions. The CD nozzle is also needed at the exit for accelerating the fuel to supersonic velocities for the MHD generator. The arrangement of the balanceofplant is shown in Fig. 3 . 1 . Thermodynamic and Heat Transfer Modelling The fuel (UFeHe mixture) is assumed to behave like an ideal gas with constant specific heats (for the purpose of this study, changes in the chemical composition of the fuel are neglected). A lumped parameter approach is selected for the thermod)oiamic and heat transfer treatment of the system. For example, it is assumed that an average value of the fuel temperature may be used to describe the heat transfer properties of the fuel inside the core at any given time. The same is true for thermodynamic properties like pressure, enthalpy, etc.. An average value is also assumed for the temperature and the heat capacity of the moderator.
PAGE 56
48 So, for example, we can write differential equations governing the dynamic behavior of the average core temperature from a simple heat balance of the form Core Heat Accumulation Rate = Core Heat Deposition Rate Core Heat Removal Rate. (31) The rate of heat accumulation in the fuel is given by MHt) C^ ^ dt where Mf(t) = instantaneous fuel (UFe+He) mass inside the core (Kg), Cpf = average fuel specific heat at constant pressure (J/KgK), and Tf (t) = average instantaneous core fuel temperature (K). The rate of heat deposition within the core is taken to be a constant fraction of the instantaneous core fission power, i.e., Qcore = 0.9 * N (32) where Qcore = heat deposition rate within the core (watts), and N = instantaneous core fission power (watts). So, it is assumed that ninety percent of the total fission energy generated within the core gets deposited within the core. The other ten percent escapes from the core in the form of prompt and delayed gammas plus the kinetic energy of the fission neutrons. Almost all of this energy is assumed to be deposited in the moderator/reflector.
PAGE 57
49 Three mechanisms of heat loss from the core are considered: radiative transfer to the wall/moderator, convective heat transfer, and the loss due to the energy transported by the flowing fuel into and out of the core. The net radiant heat transfer rate between the core and the wall /moderator Â—which are taken to be at uniform temperatureÂ— can be wrritten as (see Ref. 22 for details) 9rad = Ac O (Â£ Tf^ OC (33) where Qrad = net radiant heat transfer rate between the core and the wall (watts), Ac = surface area of the wall/beryllium moderator which is exposed to the core (fuel) (m2), o = StefanBoltzmann's constant (watts / m2K4 ), Tf= average fuel temperature inside the core (K), Tw = average wall/moderator temperature (K), 8 = average fuel emissivity, and a = average fuel absorptivity. (Implicit in the use of the above equation is the assumption that the two bodies are at uniform temperature, i.e., radiation from each body is characterized by a single average temperature.)
PAGE 58
50 If a further simplifying assumption of a = Â£ (i.e., Kirchhoffs law) is made, the above equation becomes qrad = A o e (Tf4 Tw4). (34) (Note that in actual practice, radiative heat transfer from a gas presents a far more difficult problem than is indicated by the above equation. It is complicated by the spectral and directional effects involved with gaseous absorption and emission, and also by the dependence of both a and Â£ on temperature. But, for the present work, the simplified treatment given above is presumed to be sufficient.) The film coefficient of heat transfer, H, for incompressible pipe flow is generally based on the difference between the wall temperature and mean stream temperature. For compressible flow at high speeds, it is more appropriate to base H on the difference between the actual wall temperature T^, and the 'adiabatic* wall temperature. Taw (see Ref.23 for details). Accordingly, the total convective heat transfer rate between the core and the wall/reflector can be written as Qconv = H A ( Tw Taw ) (35) where qconv = net rate of convective heat transfer (watts).
PAGE 59
51 H = effective heat transfer coefficient from fuel to the wall (watts/m2 K), and A = surface area of the wall/beryllium moderator which is exposed to the flowing fuel in the core (m2). The adiabatic wall temperature is defined by 1 +R^^M2 T 2 where T = average static temperature of the stream (K), R = recovery factor (value between 0.87 and 0.91 for subsonic pipe flow), M = Mach number of the flow, and Y = ratio of specific heats for the flowing gas (fuel). The film coefficient of heat transfer, H, is defined for turbulent flow (the flow through the core will be highly turbulent even for relatively small fuel velocities of the order of 20 m/sec) with the help of the dimensionless Nusselt number, Reynolds number and Prandtl number as follows: Nu = 0.0364 (Re * Pr)0 75 (36) where Nu = Nusselt number. Re = Reynolds number, and Pr = Prandtl number.
PAGE 60
52 (See Ref. 23 for details of the above treatment of convective heat transfer.) Neglecting the changes in the kinetic and potential energies of the flowing fluid across the core, the energy loss due to the energy transported by the flowing fluid into and out of the core can be written as fltrans Â“ r^out Cp^Tout " n^ln Cp^Tin (37) where niout(t) = instantaneous fuel mass flow rate out of the core (Kg/ sec). min(t) = instantaneous fuel mass flow rate into the core (Kg/ sec). Tin = Average fuel temperature at the core inlet (K), and Tout = Average fuel temperature at the core exit (K). This includes the energy loss associated with the thermal energy changes across the core plus the flow work changes. Equations 3.2. 3.4, 3.5, and 3.7 can be substituted in Eq. 3.1 to get a final expression for the average core fuel temperature. A similar treatment is also performed for finding the instantaneous average beryllium moderator temperature, i.e.. Be Heat accumulation rate = Be Heat deposition rate Be Heat removal rate (3 8 )
PAGE 61
53 The rate of heat accumulation in the beryllium moderator is given by Mee CpBe (dTse / dt) where Mbc = mass of beryllium moderator/ reflector (Kg), CpBe = average specific heat of beryllium (J/KgK), and Tbc = average beryllium temperature (K). The rate of heat deposition in beryllium has three components to it. The first one is the heat deposition due to the neutrons and gammas escaping from the core. This is taken to be a constant fraction of the instantaneous core fission power. So, QBe = o.l*N (39) where QBe = rate of heat deposition due to neutrons and gammas in the beryllium (watts), and N = instantaneous core fission power (watts). (That is, here it is assumed that the ten percent energy lost from the core due to escaping neutrons and gammas is given up fully in the moderator. See also Equation 3.2.) The second source of heat deposition in the beryllium moderator is the radiative transfer from the core. This part is given by Equation 3.4 above. The third source of heat transfer to the beryllium is the convective transfer from the flowing fuel inside the core to the beryllium wall. This contribution is given by Equation 3.5.
PAGE 62
54 External Loop Modelling An ideal Brayton cycle is taken to be the thermodynamic cycle. The TS diagram for the cycle is shown in Fig. 3.4. The heat addition process within the core, and the heat rejection to the heat exchanger/radiator are assumed to take place at constant pressure. The MHD generator is modelled as a conventional turbine. The performance of the turbine is then dictated by selecting an enthalpy extraction ratio and an isentropic efficiency for the turbine (these are specified as constants to the program). The values selected for the enthalpy extraction ratio and the isentropic efficiency are based on some earlier analysis of the MHD generator performed by Gerard Welch at the U of F. The enthalpy extraction ratio, together with the isentropic efficiency, can be used to obtain values for the stagnation pressure and temperature change across the turbine (MHD generator). The performance of the compressor is handled by specifying an isentropic efficiencyÂ— which is also taken to be a constant. In general, the analysis of the heat exchanger/radiator is very involved and, in the absence of detailed prior analysis, it cannot be treated in a simplistic way as is done for the MHD generator or the compressor (see discussion above). That is, the behavior of the heat exchanger/ radiator cannot be described by one or two performance indices as is done for the compressor or the MHD generator. To place the treatment of the heat exchanger/radiator at the same level as the analyses of the MHD generator and the compressor.
PAGE 63
55 o O ji Figure 34. TS Diagram for the BPGCR Power System
PAGE 64
56 and also to reduce the computational time for the transition to burst power analysis, one of two alternate simplifying assumptions can be made: either it can be assumed that the heat exchanger/radiator combination rejects heat at a constant rate or it can be assumed that it rejects a constant fraction of the heat generated within the system (the total heat includes the heat generated within the core due to fissions occurring there plus the work done by the compressor on the system). For a gas thermodynamic cycle like the Brayton cycle, the work done by the compressor is quite significant (see discussion below). Note that selection of either assumption implies that the radiator has to be of variable size (area). This is because of the fact that during the transition to burst power, the heat generated by the system (and, hence, the temperature of the working fluid) is continuously increasing, and if it assumed that the heat rejection rate is constant or is a constant fraction of the heat generated within the system, then the radiator size must be changing continuously. Such variable area radiator concepts have been proposed. For the purpose of this analysis, the second alternative is selected since it predicts an increasing heat rejection rate during the transition and, hence, is closer to the actual situation than the first assumption (i.e., in reality as the system power increases, the amount of heat rejected is also expected to increase, and dramatically so because of the T4 dependence of the radiative heat transfer).
PAGE 65
57 Bravton Cycle vs Rankine Cvrle For the purpose of this dissertation, the thermod5mamic cycle is taken to be an ideal Brayton cycle. The primary motivation for emplo)dng a gaseous cycle like the Brayton cycle was that the fuel/working fluid considered initially was a UFeHe mixture, and at typical operating temperatures and pressures existing within the system, this mixture always remained in a gaseous form. Also, initially, the possibility of a pulsed system with acoustic amplification was considered for the burst power mode of operation. It was conjectured that under these conditions, a Bra)d:on cycle might prove to be satisfactory. INSPI has also looked into the possibility of using a fuel/working fluid mixture consisting of UF4 and helium. The UF4 can be condensed out of the mixture at the temperatures and pressures existing within the system and, hence, the possibility of a hybrid Rankine /Brayton thermod)mamic cycle was also considered. Currently, INSPl is investigating fuel/working fluid mixtures of UF4 and metal fluorides in which the entire fuel/working fluid mixture is condensible under expected operating conditions. Such a fuel/working fluid mixture means that the system can operate on a true Rankine cycle. The main drawback of the Brayton cycle is the compressor work needed to compress the gas. This is significantly higher for the Brayton cycle than for a Rankine cycle. In the latter, since only an almost incompressible liquid is pumped, the typical hackwork ratio
PAGE 66
58 (which is defined as the ratio of the pump or compressor work to the turbine work) is only a few percent or less compared to 5080% for a Brayton cycle (for an ideal Bra)d;on cycle it is of the order of 50% but in actual systems, it can be as high as 80%). A second drawback of the Brayton cycle, which is very significant for space power systems, is the temperature at which heat rejection occurs. In a Rankine cycle, the heat rejection can take place at constant temperature (condensing vapor) whereas in a Brayton cycle it occurs at a varying temperature. If we were to determine an "effective" constant temperature for heat rejection for a Brayton cycle, it turns out to be much lower than the constant temperature at which the Rankine cycle rejects heat (provided we place comparable thermal limits on the Brayton and Rankine cycles). This factor can dramatically increase the size of the radiator. If the frictional losses in the system are neglected, and if the efficiencies of a Brayton cycle and a Rankine cycle are compared (for roughly the same operating conditions), the Rankine cycle can actually turn out to be poorer than a Brayton cycle (2226% for a Brayton cycle compared to roughly 22% for a Rankine cycle) (numbers are based on some earlier work done by Gerard Welch at the U of F)). However, we cannot neglect frictional loses, and making up for these in a Rankine cycle costs almost nothing; the cost in a Brayton cycle can be very high because of the high hackwork ratio. With the frictional loses included, the Brayton cycle efficiency turns out to be lower than that for the
PAGE 67
59 Rankine cycle. This lower efficiency also translates into a bigger radiator mass (since a lower efficiency implies more heat rejection for the same electrical power output). Fluid Flow Modelling As mentioned earlier, the fuel (UFeHe mixture) is assumed to behave like an ideal gas with constant specific heats. The fuel is also assumed to maintain its integrity throughout the transition to burst power (i.e., no dissociation of the UFÂ© is taken into account). Even after the commencement of the transition, it takes some time before the core temperature and pressure builds up to levels where significant dissociation of UFe can occur. Significant UFe dissociation occurs only at temperatures above 1600 K at atmospheric pressure. With increasing gas pressure, dissociation becomes significant only at higher temperatures. For e.g., with a UFe partial pressure of about 10 atm., the fuel can exist without much UFe dissociation up to about 2200 K. The fuel temperature profile during the transition is given in Fig. 3.21 later in this Chapter. From this figure, it can be seen that, at least for the first 50 sec. of the total 100 sec. transition time, the fuel temperature remains below the 2200 K value. (UFe may be able to exist without significant dissociation at much higher temperaturesup to about 3000 K, provided fluorinating agents like CIF 3 or BiFa are added to the fuel (Ref 24)).
PAGE 68
60 If the fuel is in the form of UF 4 . the dissociation is not at all a concern even at temperatures as high as 30004000 K. Neutronically, it doesn't matter whether the fuel is in the form of UF 4 or UFe. The extra fluorines have no significant neutronic effect. So, the neutronic analysis done here remains valid whether UFe or UF 4 is used as the fuel. Core hydrodynamics is represented by fully developed, onedimensional, turbulent flow (as mentioned earlier, the fuel flow through the core will be highly turbulent even for relatively low fuel velocities of the order of 2050 m/sec). Frictional effects are neglected. Isentropic flow is assumed for the flow through the inlet and exit nozzles so that stagnation temperatures and pressures of the fuel will adequately describe the flow of the fluid through the nozzles. For example, the maximum fuel flow rate through the exit nozzle is given by the relation u, _ 0.6847 Po A i^max TTo R (To)^/2 (310) where nimax = maximum fuel mass flow rate (Kg/sec), A = exit nozzle throat area (m 2 ) , Po = core stagnation pressure (Pa), and To = core stagnation temperature (K).
PAGE 69
61 (This maximum flow occurs under choked flow conditions, i.e., when the Mach number at the throat is unity). The exhaust nozzle is assumed to be under choked flow conditions always (flow through the convergingdiverging nozzle at the core exit has to be choked to obtain supersonic flow going into the MHD generator). The converging nozzle at the core inlet is operated under nonchoked conditions for reasons given earlier in the geometry modelling section. Neutron Kinetics Modelling It was mentioned earlier that it is assumed that negligible decomposition of UFe takes place during the transition to burst power. Apart from this, it is also assumed that the fuel composition remains constant while going through the external loop. This implies that the UFe mole fraction of about 93% and the 85% uranium enrichment remains the same. It is also assumed that there is no significant buildup of strongly neutron absorbing fission products. These assumptions should not result in too much error since the fuel bumup and fission product buildup are going to be very small during the short transition time. To model exactly the time behavior of the neutron population within a reactor, use of the timedependent Boltzmann's neutron transport equation has to be made. This is all the more true for the type of gas core reactor that is anal}^ed in this dissertation because of the fact that within the core proper, almost no neutron moderation
PAGE 70
62 (scattering) is taking place, and the core is only of the order of a mean free path or less in size. The only significant interaction taking place within the core is the absorption of thermal neutrons returning to the core from the reflector/ moderator. But, from a practical standpoint, implementing this model can be very tedious and computationally intensive, and is attempted only in major time dependent transport theory computer codes. Some previous studies conducted at the U of F indicate that, for gas core reactor systems which posses a significant degree of grayness, multigroup diffusion theory is generally adequate for static calculations provided good fast and thermal group constants are used (Ref. 25) (grayness for an externally moderated gas core reactor is defined as the ratio of the thermal neutron current into the core from the moderating/reflecting region to the value of the thermal neutron flux at the coremoderator/reflector boundary). For small transients (i.e., small reactivity variations), we can assume that timedependent, multigroup diffusion theory is adequate to describe the dynamic behavior of the gas core reactor. Sometimes, even this model is too detailed for practical implementation in reactor kinetics analysis due to excessive computational requirements. For this reason, in this dissertation (as is the common practice for many preliminary level reactor kinetics studies), it is assumed that the neutron dynamic characteristics of the gas core reactor can be adequately described by the average neutron density, 'n', i.e., the
PAGE 71
63 spatial dependence of the neutron population is neglected. Often, in neutron kinetics studies, 'n' is assumed to represent only the neutrons in the energy range where most of the fissions occurÂ— for example, thermal neutrons in a thermal reactor. This type of analysis is by no means exact but, for preliminary scoping purposes, this is adequate. It is capable, for example, of providing the time dependence of the gross reactor power level during transients. With the above assumption, the d 5 mamic behavior of the neutron density within the core can be represented by dn(t) dt p(t) (3eff r n(t) + X i where n(t) = core neutron density at time t (#/cm3), p(t) = time dependent system reactivity. Peff = effective delayed neutron fraction. (3.11) 1* = mean neutron generation time (sec), Ci(t) = effective delayed neutron precursor density for the ith delayed neutron group (#/cm3), and A .1 = decay constant for the ith delayed neutron precursor group (seci). To complete the description of the neutron dynamics of the reactor, an equation for the ith delayed neutron precursor density, Ci, also has to be provided. For the gas core reactor analyzed here, this
PAGE 72
64 precursor density equation is quite different from that for a conventional, solid core reactor. Since this is a circulating fuel reactor, not all of the delayed neutron precursors will decay within the core and thereby contribute to the total core neutron population. Some of the precursors are going to be swept out of the core (along with the fuel) before they can decay. Now two things can happen to these precursors which have left the core without decaying; they can either decay outside the core (so these neutrons are lost forever, and will not contribute to the neutron population within the core) or they can reenter the coreÂ— after going through the external loop without decaying. In either case, the effective fraction of delayed neutron for this type of reactor is going to be less than that for a stationary fuel reactor. The precursor density equation taking into account the fuel circulation is dCj(t) ^ Peff n(t) _ ^ , _ ci(t) ci(txi) exp (Xi ^j) dt 1* ^ ^ (312) for i = 1,2 6. where Xc = mean core residence time for the fuel (sec), and Ti = mean external loop circulation time for the fuel (sec) (All other variables are as defined in Equation 3.11.) The third term on the right hand side of the above equation accounts for the loss of the delayed neutron precursors from the core due to
PAGE 73
65 flow, and the last term on the right hand side takes care of the reentry of the precursors at the core inlet after going through the external loop. So, the effect of the fuel flow outside the core is to make the reactivity needed (and, hence, the fuel mass needed inside the core) for steady state operation larger than the value needed if it were a nonflow system. In Equation 3.11, 'n' represents the instantaneous neutron density within the core. But it can represent any integral or volume averaged property of the reactor that is proportional to the instantaneous neutron density. For example, it can be redefined so that it represents the total number of neutrons within the core at any given instant or the fission rate or the total core power at any given instant, etc.. For convenience, the dependent variable in Equation 3.11 is defined to be the instantaneous power, N, of the reactor. Now the precursor density, Ci, used in Equations 3.11 and 3.12 has to be replaced by a new quantity Ci defined as follows: Ci = WfVXfCiv (313) where Cl = effective delayed neutron precursor density of the ith delayed neutron group (#/cm3), Wf = usable energy produced per fission (MeV/fission),
PAGE 74
66 Zf = average thermal macroscopic fission cross section (cm1 ) , V = total core volume (cm3), and V = neutron velocity (cm/sec). As a final simplification, instead of the usual six delayed neutron groups, only one effective group of delayed neutrons is considered, i.e., the summation term in Equation 3.11 is replaced by one term which represent an "effective" delayed neutron group. Also, Equation 3.12Â— which is a set of six equations, is now replaced by one equation for the "effective" delayed neutron group. This simplification is also done so as to reduce the computational time. The computer program developed can equally well handle all six delayed groups (with a corresponding increase in the computational time). There are different ways of defining an "effective" one delayed group of neutrons. The recipe used in this dissertation for obtaining the value for J3, and the corresponding value for X is as follows; .6 = I^i (314) and (I Pi/ Xi) i (315) for i = 1,2 6
PAGE 75
67 where = fraction of the ith delayed neutron group emitted per fission, and = decay constant for the delayed neutron precursor for the ith delayed neutron group. The neutron generation time, T, used in Equations 3.11 and 3.12 is an important factor, and is defined as the average time between the birth of a neutron and its subsequent absorption inducing fission. It is related to the prompt neutron lifetime, T of the reactor by the relation 1* = 1 / keff where 'keff is the effective neutron multiplication factor of the reactor. The prompt neutron lifetime (which is defined as the average time between the birth of a neutron and its subsequent removal from the systemÂ— either by absorption or by leakage) is a very difficult parameter to estimate anal 5 dically since it depends in a complex manner on a number of system variables like geometry, fuel enrichment, reflector thickness, etc. So, for the purpose of this work, T for the system is determined by independent Monte Carlo calculations, and the value thus obtained is used as an input to the program.
PAGE 76
68 Programming Aspects The computer program, DYNAM, developed as a part of this work for analyzing the transients, is written in standard FORTRAN 77. It can run either on a PC or a mainframe (the program calls for certain timing subroutines which are only available on IBM mainframes. These call statements need to be deleted from the program before tr 3 dng to run it on a PC or on a nonIBM mainframe). Simple finite difference techniques are used for solving most of the relevant equations like those for core temperature, beryllium temperature, core fuel mass, etc.. The coupled neutron kinetics equations are solved using a fourth order RungeKutta method (this method is selected because it is quite easy to implement, plus it has good numerical stability characteristics). The only difficult term to handle in the neutron kinetics equation is the precursor reintroduction term, i.e., the last term on the right hand side of Equation 3.12Â— which represents the precursors which are reintroduced into the core after going through the external loop. This term is complicated because, at any given instant, it depends on the precursor density which had existed within the core at a time Ti (the mean external loop residence time) before. When the core power is varying rather slowly, a reasonable approximation would to be take
PAGE 77
69 Cl (tTi) = Cl (t) This is certainly not true during the transition to burst powerÂ— where the power goes from almost zero to a few hundred megawatts in about 100 sec. So, in DYNAM, the precursor reintroduction term is handled numerically as described below. The precursor density at each time step is stored in a stack, and, at any instant, the value of 'C which existed at a time li before is used to estimate the precursor reintroduction term. Once a value C(tXi) is used during a time step, it can be discarded from the stack of stored values of 'C, That way, at any given instant, the number of values of 'C that need to be stored is only Xi / step where 'step' is the (constant) time step size used in the numerical computation. Sometimes this can create memory limitation problems while running on a PC. In this case, the precursor density can be averaged over a fixed number of time steps before storing it. The number of times steps over which the precursor density is averaged depends on the transient. If it is a very fast transient, one might want to average the density over a small number of time steps, whereas, for slower transients, the density can be averaged over a large number of time steps without incurring serious error.
PAGE 78
70 Testing of the Computer Program Before the developed computer program, DYNAM, can be used for analyzing the transition to burst power, it has to be tested or benchmarked. The complete programwhich includes the neutron kinetics, the thermodynamic, the heat transfer, and the fluid flow aspects of the gas core reactor, cannot be benchmarked against any standard codes since there are no existing codes which are capable of analyzing transients for such a novel type of reactor (novelty includes the combination of gaseous fuel in a closed circulating loop, highly compressible fuel which also serves as the working fluid, etc.). So, it was decided to test individual sections of DYNAM against some standard codes or against analytical resultsÂ— where possible. One of the major subroutines in DYNAM is the one which estimates the effective neutron multiplication factor of the system (or the system reactivity). This subroutine accepts the reactor dimensions (core radius and beryllium thickness) and the fuel atom densities as input, and it estimates the system keff. The system keff needs to be known at every instant during the transition to burst power. Since the system keir has to be estimated a large number of times during the transition analysis, this subroutine was written giving speed of calculation priority over accuracy. The subroutine mocks up the cylindrical core with a spherical geometry, and then performs a tworegion (core and moderator/reflector), twogroup (fast and
PAGE 79
71 thermal) diffusion theory kgff calculation. The microscopic group constants needed for the keff calculation are provided as input to the subroutine. These group constants are obtained from independent, twogroup transport theory calculations using XSDRNPM. These microscopic constants are found to be not very sensitive to small changes in fuel atom densities within the core or to changes in the core dimensions. But, if either of the two (geometry or fuel atom densities) is changed significantly, the group constants need to be reevaluated using XSDRNPM. The effects of temperature on the microscopic constants are handled via feedback coefficients. These coefficients are given in Tables 3.6 and 3.7. The system keir values calculated by the DYNAM subroutine are compared against those estimated by XSDRNPM (for the same geometry and fuel atom densities). Comparisons are made against twogroup diffusion theory and twogroup transport theory results from XSDFdMPM. The results of the comparisons are shown in Table 3.9. The comparisons are done for various fuel atom densities within the core. As can be seen from the table, in all cases, the error involved is less than a percent and, hence, is quite acceptableÂ— especially considering the fact that the subroutine is quite fast. It takes less than 0,01 sec. to estimate the keff of the system to a convergence level of better than 104 on an IBM 3090 mainframe. Compared to this, a
PAGE 80
72 Table 39. Comparison of keff Values Obtained Using the kecr Subroutine of DYNAM and the XSDRNPM Code for Different Fuel Density Factors. 1 Core Fuel Density Factor keff Percentage Difference Between DYNAM keff and 2 DYNAM 3 2Group XSDRNPM (Diffusion Theory Option) 4 2 Group XSDRNPM (S4P^ Two Group XSDRNPM (Diff. Theory) Two Group XSDRNPM [S^P^ 0.2 0.6990 0.7025 0.7000 0.498 0.143 0.5 1.0356 1.0411 1.0384 0.528 0.269 0.8 1.1760 1.1789 1.1758 0.246 +0.017 1.0 1.2310 1.2337 1.2302 0.219 +0.065 1 . A fuel density factor of unity corresponds to the following atom densities (in atoms/barncm); U235 = 1.458 E5 U238 = 2.548 E 6 F = 1.027 E4 He = 2.275 E4 2. Two group, two region, diffusion theory calculation on a spherical core with 72. 1 cm. core radius and 50 cm. beryllium moderator/reflector thickness. 3. Two group, two region diffusion theory calculation on a spherical geometry of the same dimensions as given above. 4. Two group, two region transport theory calculation (S 4 P 3 ) on a spherical system of the same dimensions as given above.
PAGE 81
73 twogroup, tworegion XSDRNPM keff calculation on a similar system takes about 0.5 sec. The fact that the DYNAM results are closer to the two group transport theory results is purely coincidental. In general, it was found that diffusion theory predicts a keir which is higher than the keff prediction of transport theory. So, the XSDRNPM two group transport theory results are lower than the XSDRNPM diffusion theory results. It so happens that DYNAM predicted keff values are also slightly lower than the XSDRNPM two group diffusion theory results and, thus, actually a little closer to the two group transport theory results. To complete the benchmarking of the keff subroutine of DYNAM, the fast and thermal group fluxes obtained using the keff subroutine are also compared against those obtained using two group transport theory calculations (XSDRNPM). The results of this comparison are given in Figs. 3,5 and 3,6. As can be seen from the figures, except for the case of fast flux within the core, DYNAM predictions of the thermal and fast fluxes are in agreement with the predictions of XSDRNPM. DYNAM underpredicts the fast flux within the core by roughly 10%. The fluxes shown in these figures (group 1 and group 2) are normalized such that the total fission neutron source in the core is one neutron per second (it should be mentioned here that the agreement between flux shapes predicted by the diffusion theory and the transport theory is not always this good. For certain fuel density regimes, the agreement
PAGE 82
2.00E4 74 S9xn[j (;sBj) I dnojo p9zixBUU0N Figure 35. Comparison of Radial Flux (Fast) Profile Within the BPGCR
PAGE 83
2.00E4 s ll
PAGE 84
76 between Sn and diffusion theory fluxes can be very poor even though the keff values are in reasonable agreement). A second major subroutine in DYNAM is the one which solves the coupled point reactor kinetics equations. This subroutine uses a fourth order RungeKutta method to solve the coupled neutron and precursor density equations. Two separate methods were used to check the validity of this subroutine. First, the effectiveness of this subroutine in predicting the neutron density behavior after step insertions of reactivity was compared against analytical results. For two different step reactivity insertions, the results obtained using the developed subroutine are compared against the "exact" analytical results (see Ref. 26), and also against the 'prompt jumpÂ’ approximation results (the prompt jump approximation assumes that immediately following a step change in reactivity, the delayed neutron precursor concentration remains unchanged over the time of the sudden rise or drop in flux). The "exact" analytical result is for one delayed neutron group: it gives good results when the reactivity step is no greater than about 0.3 J3 (see Ref.26 for the limitations of this "exact" anal}d;ical expression). Figures 3.7 and 3.8 gives the results of this comparison. Again, it can be seen that the results from the developed subroutine agree quite well with analytical predictions. Next the developed subroutine is compared against the standard point reactor kinetics code ANCON (Ref.27). For the purpose of
PAGE 85
116 77 o 9i cn M M o cd V Qj V % V a H jamoj jo^3bbh Figure 37. Benchmarking of DYNAM Against Analytical Methods Positive Reactivity Insertion Case
PAGE 86
p(t) = '0.00065 p(0) = 0.0 N (0) = 100.0 watts 78 JdAtOd JO^OBBH Figure 38. Benchmarking of DYNAM Agaist Analytical Methods Negative Reactivity Insertion Case
PAGE 87
79 comparison, three different ramp reactivity insertion rates were selected. Either six delayed neutron groups or one delayed neutron group treatment can be used in ANCON. Figures 3.9 and 3.10 compare the results from DYNAM against six delayed neutron group ANCON results (note that the one delayed group ANCON results agree exactly with the DYNAM results). Since ANCON cannot deal with circulating fuel type reactors (it assumes that the fuel is stationary), for the purpose of comparing the results of DYNAM with those of ANCON, the fuel circulation effects (for example, the loss of delayed neutron precursors from the core, the delayed neutron precursor reintroduction, etc.) are zeroed in DYNAM. Results The developed computer program, DYNAM, is now used to study the transition to burst power. As mentioned before, the program is very versatile and can be used in its present form or with very little modification to study different aspects of the rapid transition. In the following sections, a few selected capabilities and applications of the program are demonstrated. First, two separate reactivity insertions are applied to the core to bring it up from almost zero power to full power in roughly 100 seconds. During this transition, crucial system parameters like core pressure, core temperature, etc., are estimated by DYNAM. Based on the behavior of these parameters, system design can be developed or modified. For example, if one finds that, for some desired operating
PAGE 88
Relative Neutron Density 80 Figure 39. Benchmarking of DYNAM Against ANCON Case 1
PAGE 89
Relative Neutron Density 81 Time After Ramp Reactivity Insertion (sec) Figure 310. Benchmarking of DYNAM Against ANCON Case 2
PAGE 90
82 power level, the average core temperature is higher than some acceptable value, the system design can be modified, e.g,, to increase the fuel flow rate through the core so as to reduce the average core fuel temperature. This increased fuel flow rate through the core might require changing the ratings of the pumps or even changing the core inlet and exit nozzle areas. The behavior of these crucial system parameters can also be used in selecting the best possible reactivity insertion pattern that is capable of taking the system from zero power to full power in a safe manner. The way DYNAM is set up currently, to bring about the transition to burst power, one has to specify the desired reactivity variation (this can be almost any conceivable reactivity variation). The program Avill then estimate the fuel mass variation needed within the core to bring about this reactivity variation. It also determines the fuel flow rates needed in order to obtain the necessary fuel mass inside the core, and the core inlet pressure is automatically adjusted to give the required fuel flow rates. There are many types of reactivity variations that can be applied to take the core from zero power to full power. For the purpose of demonstration, two different reactivity insertions are selected. The first one is a pure linear (or ramp) reactivity variation. In fact, as shown in Fig. 3.11, it is a combination of two linear reactivity variations. To simulate the conditions in an actual transition scenario, the core is initially assumed to be far subcritical (ptot = 1.0). This
PAGE 91
Figure Not to Scale 83 .. 8 o o (N o o a o ed & tn o % o rt/w ) iC;iAnoB3H Hcx)
PAGE 92
84 would then correspond to a situation where the fuel is just being introduced into a previously unfuelled (voided) central chamber. During the first part of the reactivity variationÂ— when the core is subcritical, a very steep reactivity variation of 0.02 Ak/k per second is introduced into the core (for a U235 fueled system with a delayed neutron fraction of about 0.0065, this number translates to roughly $3.0 per second). Once the system reactivity becomes zero (keff = 1.0), the slope of the linear reactivity variation is reduced to 0.00015 Ak/k per second (roughly $0.02 per second). Figures 3.12 through 3.16 show the variation of some of the important system parameters during the rapid transition to about 100 MWt as computed by DYNAM. The core power is set at an initial value of 1000 watts within the program. This might correspond to a situation where the bimodal system is getting ready to go to a high power mode of operation from a low power mode, and some power is generated within the BPGCR due to neutrons reaching it from the surrounding PGCRs. The behavior of the core power (Fig. 3.12) can be explained based on the behavior of the reactivity variation alone as follows. Initially, when the core is far subcritical, the core power falls as time increases even though the system reactivity is increasing. This happens due to the fact that the system is so far subcritical that it is not able to maintain a steady state neutron population. Then the total system reactivity reaches a point (still subcritical) where the rate of fuel mass increase within the core forces the power to go up (note
PAGE 93
85 that the core power is given by an expression which involves a product of the flux and the fuel atom densities within the core. So, even though the flux is still falling (since the core is still subcritical), the rate of increase of the fuel atom densities compensates for this effect and the net result is a power that increases with time). Then at 50 seconds after the initiation of the transition, the system reactivity becomes exactly zero (critical stationary system), and, at this point, within DYNAM, the positive reactivity insertion rate is reduced to a much lower value (see Fig. 3.11). The system power now drops slightly due to fact that, for this circulating fuel system, there is some loss of delayed neutrons while going through the external loop and, even though the stationary system reactivity is zero at time t = 50 seconds, the circulating system reactivity is still negative (system multiplication factor is less than one). Eventually, at t = 85 seconds, the circulating system becomes critical and further reactivity addition leads to a rapid power increase. To understand the behavior of other system parameters (like core temperature, core pressure, etc.) the effect of the balance ofplant also has to be taken into account. The effect of reactivity variation alone cannot explain the behavior of these system variables. For example. Fig, 3.14 shows that the core temperature increases from the very beginning of the transition even though the core power is dropping during the initial phases of the transition. The
PAGE 94
lE+9= 86 o jBuuaqx aiOQ Figure 312. BPGCR Thermal Power Vs. Time Case I
PAGE 95
OOT 87 o (3s) 9100 axi:Â» ssbjv pnj I^^ox Figure 313. BPGCR Total Fuel (UF6+He) Mass Variation With Time Case I
PAGE 96
1600 88 (H) 3Jn:^BJ3duiax dzoQ sSvjqav Figure 314. BPGCR Average Core Temperature Variation with Time Case I
PAGE 97
30.0 89 (UI^B) 3JnSS3J
PAGE 98
160 90 o u V Â« c 0 (tf a p, p 1 Â•M CÂ» ll Â« V B Â•pN H (s/3h) Momtio ssBji\[ xaiva aioo Figure 316. BPGCR Fuel Mass Outflow Rate Variation with Time Case I
PAGE 99
91 reason for this behavior of core temperature can be directly traced back to the compressor work done on the fuel/ working fluid. As mentioned earlier, since this is a Brayton cycle, the compressor work needed to compress the gaseous fuel is considerable. This effect dominates during the initial phase of the transition and the core power (fission power) has very little effect on the core temperature. The work done by the compressor depends on the fuel mass flow rate and, hence, initially, during the rapid insertion of reactivity (addition of fuel to the core), the work done by the compressor is also increasing. Only towards the end of the transition, when the core power becomes significant (above a few megawatts), does the core power starts influencing the behavior of core parameters like core pressure, core temperature, etc.. From Fig. 3.16, we can see that most of the fuel mass flow rate increase occurs during the first 50 seconds (except for the sharp increase towards the end) of the transition. This behavior of the fuel mass flow rate gets reflected in the fuel temperature profileÂ— where the core temperature pretty much comes to a steady state after the first 50 seconds or so. The discontinuity in the core parameters around 50 seconds after the initiation of the transition can be attributed to the discontinuity in the reactivity insertion rate at this point in time (i.e., a change in the reactivity insertion rate results in a change in the fuel mass inflow rate and, hence, a change in the compressor work output). The core temperature profile does not show any sharp discontinuity at 50
PAGE 100
92 seconds because of the damping that occurs in the external loop (due to the heat exchanger/radiator). The type of reactivity insertion described above (pure ramp insertion) is not a good one to employ during a real transition. As can be seen from Figures 3.12 through 3.16, most of the core power increase (and, hence, the increases in the core temperature, core pressure, etc.) occurs very rapidly during the last few seconds of the transition. This means that some type of drastic controlling action needs to be taken at the end of the transition period to prevent power overshoot. Also, most of the time, the safety features engineered into a reactor are designed to prevent such rapid power increases, and will cause the reactor to be scrammed. So. instead of using just a pure ramp insertion of reactivity, a specific reactivity variation, as shown in Fig. 3.17, was next selected for the purpose of demonstration. This reactivity variation is designed to produce a fast power increase to a predetermined constant power level. As shown in Fig. 3.17, the reactivity variation is a combination of two linear ramps (with different slopes) and an exponentially decaying reactivity variation. By adjusting the peak value of the core reactivity (Pmax) and the time at which the exponentially decaying reactivity variation begins (ti), the steady state power level attained can be varied (see Ref. 28 for details of this type of reactivity variation. Ref. 28 also has details about other types of
PAGE 101
0.0074 93 (3I/3IV ) HOOda Figure 3.17. BPGCR Input Reactivity Variation With Time Case II
PAGE 102
94 reactivity variations that are capable of producing fast monotonic power increase to predetermined constant power levels). Figures 3. 18 through 3.23 show the variation of the most important system parameters during the rapid transition as predicted by DYNAM. Again, the behavior of the system parameters can be explained as a combined effect of the type of reactivity variation introduced into the core during the transition, plus that of the balanceofplant, as was done for the previous case. Another effectÂ— which has some significant effect on the dynamics of the systemÂ— can be noted from the transition scenario described above. The beryllium temperature hardly changes during the entire transition. This is the case even when there is no heat removal mechanism provided for the beryllium moderator. In the program, the beryllium temperature is initialized as 800 K. This might represent a case where the beryllium temperature is maintained at 800 K during the low power mode of operation (see discussion in Chapter 1). The sluggish behavior of the beryllium moderator temperature is due to the large thermal inertia of the moderator. One outcome of this sluggish behavior of the beryllium temperature is that the reactivity feedback, due to the moderator temperature change, can be safely neglected during the rapid transition. It needs to be taken into account only for dynamic analysis of prolonged operational transients at high power levels.
PAGE 103
lE+9 J3MOJ IBUUailX 3IOO
PAGE 104
96 (35E) 3JOO am ssbjvi lanj ib^ox Figure 319. BPGCR Fuel (UF6+He) Mass Variation with Time Case II
PAGE 105
97 u 4 ) a o (Q a S CO Â«> % V a (ui^B) djnssajj ajoo a^BJdAv Figure 320. BPGCR Average Core Pressure Variation with Time Case II
PAGE 106
2200 98 (H) ajn^Bjadraax 3JOO b^bjbav Figure 32 1 . BPGCR Average Core Temperature Variation with Time Case II
PAGE 107
180 99 o V 0) a o 0< :3 I M (/) U 4 ) 4 ) B H (s/^H) Momno ssbi\i lan^ ajoo Figure 322. BPGCR Fuel (UF6+He) Mass Outflow Rate Vs. Time Case II
PAGE 108
1000 100 o u V (Â» C o Â•M CS a A M (Â» u % o a H (H) ajn^Bjadtaax io^Bjapoi\i umni^aa Figure 323. BPGCR Average Moderator Temperature Vs. Time Case II
PAGE 109
101 The predictions of DYNAM during the transition can be very helpful in designing the system or modifying the design of the system. An example of this capability of DYNAM can be seen by considering Figure 3.21. Here the variation of the average core temperature is given when the reactivity insertion shown in Fig. 3.17 is applied to the system. Now suppose it is decided that the steady state (at the end of the rapid transition) mean core temperature is not high enough. One way of increasing this average temperature (without increasing the core power) is to decrease the fuel mass flow rate through the core. This effect is shown in Fig. 3.24Â— where the core average temperature variation is given for the same reactivity insertion shown in Fig. 3.17, but, with a core exit throat area of 0.03 m2. For the cases shown in Figures 3.18 through 3.23, the throat area of the convergingdiverging nozzle at the core exit is fixed at 0.05 m2. This area, to some extent, dictates the fuel mass flow rate through the core (see Eq.3.10). (Note that at the core exit, the flow is kept choked so as to obtain supersonic velocities at the MHD duct inlet). It can be seen that the steady state core temperature has now increased to about 2400 K from the previous value of about 2000 K. Figure 3.25 shows the new fuel mass outflow rate from the core. Comparing Fig. 3.25 with Fig. 3,22, it can be seen that the fuel mass flow rate is now lower (resulting in a higher core temperature).
PAGE 110
2500 102 o V 0} C o Â•*4 M cd S (0 M V 4) a (H) 9ni:(BJ9diii9x aJOQ s^bjbav Figure 324. BPGCR Average Core Temperature Vs. Time
PAGE 111
120 103 u V m (3 O *4 M eS a A s M CO Ki V % V a H (s/ 3 h) Mou;no ssbj^ lan^ aJOQ Figure 325. BPGCR Fuel (UF6+He) Mass Outflow Rate Vs. Time
PAGE 112
CHAPTER 4 BPGCR STABILITY ANALYSIS Introduction Investigation of the inherent stability constitutes an essential part of a reactor design process. This is particularly true for a new type of reactor system such as the BPGCR concept being investigated in this study. A fundamental consideration in such an investigation is to determine whether the system possesses inherent selfdestruction tendencies, i.e., to determine whether a given djmamical system (in this case, the reactor) will eventually return to an equilibrium state after a finite perturbation. A reactor is clearly unstable if a reactivity perturbation can cause an unlimited increase in neutron flux or power output. This is true whether the increase in power or neutron flux occurs continuously or as a series of oscillations of increasing amplitude. For the system to be considered inherently stable (i.e., without the help of any external controlling action), a perturbation during steady state operation must eventually lead to the same or another equilibrium state of the system. The transient behavior in going from one steady state to another should also be satisfactory, i.e., the power output (or flux) should not undergo a series of oscillations of unacceptable magnitude before settling down. 104
PAGE 113
105 Several approaches may be used for stability analysis of a reactor system. A complete study of power reactor dynamics would take into account the inherent nonlinearity of the reactivity feedback. This nonlinearity is due to the dependence of the system reactivity on core power. The transient response of a nonlinear system can be obtained by solving simultaneously the system equations. As an example of this procedure, using the computer program, DYNAM (described in Chapter 3), the behavior of the BPGCR after a reactivity perturbation can be studied. DYNAM numerically solves the nonlinear neutron kinetics, thermodynamics, and heat transfer equations of the system. The power output versus time, as predicted by DYNAM, for a single chamber BPGCR system, initially operating at a steady state power of roughly 100 MWt, following a one dollar ($1) external positive step reactivity insertion is shown in Fig. 4. 1 . Note that for a circulating fuel reactor, the reactivity which corresponds to a dollar is not the same as that for a stationary fuel reactor. Since the effective delayed neutron fraction is reduced in a circulating fuel reactor due to the loss of delayed neutrons that are emitted external to the core, the reactivity that corresponds to a dollar is given by 'a J3' , where 'a' is called the delayed neutron attenuation factor, 'a' is a function of both core residence time and external loop circulation time of the fuel (see Eq. 4.4). The external loop circulation time is held constant at 0.1 second for this, as well as for most of the other transients considered in this chapter. (In an actual system, the external loop circulation time
PAGE 114
Initial Core Power = 100 MWt $1 Positive External Reactivity Step (pext = t0.0029) 106 jamoj ma^si^s IB^ox Figure 41 . BPGCR Thermal Power Vs. Time After a One Dollar External Positive Reactivity Insertion
PAGE 115
107 depends on the length of the external loop and the average fuel velocity while traversing the external loop. Typical values for the external loop circulation time should only be a few tenths of a second, and, for realistic values of external loop length and fuel velocities, a second is probably an upper limit for the external loop circulation time.) The core fuel residence time is of the order of 0.08 sec. The system reactivity is a combination of externally applied reactivity, plus the reactivity due to the fuel mass inside the core. No feedback effects due to beryllium moderator temperature change, or Doppler effect due to fuel temperature change are considered for this perturbation. The beryllium temperature takes a long time to change because of its large thermal inertia and, hence, the feedback effect due to this change is expected to be negligible for transients that are but a few seconds in duration. The fuel temperature (Doppler effect) coefficient of reactivity is neglected because of its small magnitude compared to the fuel mass (or density) feedback. For the illustrated transient, the inlet stagnation pressure is kept constant. During normal operation, the core inlet stagnation pressure is constantly monitored and is adjusted so as to obtain the necessary fuel mass inside the core. However, for this external step reactivity insertion, the value of the inlet stagnation pressure is frozen at its initial steady state value, and is kept constant during the transient. As shown in the figure, as soon as the external step reactivity is applied, the core power begins to rise. This results in an increase in
PAGE 116
108 the core temperature and pressure. But since the core inlet stagnation pressure is kept constant, the increase in core pressure has two effects; first the fuel mass flow rate at the core inlet starts decreasing (the core inlet is always maintained at nonchoked flow conditions) and second, at the same time the fuel mass flow rate at the core exit starts increasing. These two effects drive the total core reactivity down and the core power starts to decrease. The rate of power decrease is halted (or reduced) only when the core pressure has fallen enough to make the constantly held inlet stagnation pressure high enough (with respect to the core pressure) to make the total core reactivity (fuel mass) go up again. It can be seen from the figure that the strong fuel mass (or density) feedback has the effect of bringing the system under control rather quickly (in less than a second) even after a hefty reactivity perturbation ($1 step insertion). The amplitude of oscillation is also very modest (less than 50% of the initial steady state value) during the initial phase of the transient. The system behaves in a safe manner after the reactivity perturbation, i.e., no runaway power output ensues after the perturbation, even though no external controlling action is provided (except for maintaining the core inlet stagnation pressure constant during the transient). The BPGCR internal reactivity (reactivity due to the fuel mass inside the core) is plotted versus time for this transient in Fig. 4.2. During the entire transient, the externally applied reactivity is held constant at one dollar.
PAGE 117
0.0040 109 ( n/yiv) iBtua^ui ajoQ Figure 42. BPGCR Internal Reactivity Vs. Time After a One Dollar External Positive Reactivity Insertion
PAGE 118
no From Fig. 4.2, it can be seen that the net (internal plus external) system reactivity remains positive during the entire 2 second transient, and yet, the system power is decreasing. This is happening because of the loss of the delayed neutrons due to fuel circulation, i.e., the decay of some of the delayed neutron precursors while in the external loop effectively reduces the delayed neutrons available to the system for multiplication. So, the circulating fuel reactor reactivity has to be positive (and not just zero) for the system power to remain steady or to increase. From Fig. 4.1, it can be seen that an external positive reactivity step results in a lower steady state power compared to the steady state power level before the external reactivity insertion. It might then be inferred that a negative external reactivity step would result in a steady state power level that is higher compared to the power level before the reactivity insertion; this would be a very undesirable response. So, to study the effect of a negative reactivity insertion, a negative one dollar reactivity step is applied to the BPGCR. The resulting power behavior is given in Fig. 4.3. It can be seen that, for this case also, the final steady state power level is lower than the power level before the perturbation. The internal reactivity variation is shown in Fig. 4.4. Comparing this figure with Fig. 4.2, it can be seen that, during the initial oscillations, the behavior of the internal reactivity after a negative external step reactivity insertion is almost a mirror image of that after a positive external reactivity step insertion.
PAGE 119
1.2E+8 111 O V 0 ) a o t <Â» C M *5 o ed V H e V M EH) V 6 Figure 43. BPGCR Thermal Power Vs. Time After a One Dollar External Negative Reactivity Insertion
PAGE 120
08000 112 o 9i Â•M CO *3 e 4 > 4) 4 > a ( 31/317) iC;AI^3B3H IBUJO^UI OJOO Figure 44. BPGCR Internal Reactivity Vs. Time After a One Dollar External Negative Reactivity Insertion
PAGE 121
113 This response of the BPGCR is quite different from that of a conventional, solid core (stationary fuel) reactor, i.e,, the BPGCR steady state power cannot be increased merely by an external positive reactivity addition. Although the behavior exhibited in Figs. 4.1 and 4.2 may be a very desirable one during accident situations (abrupt insertion of positive or negative reactivity), it can pose some serious operational problems during normal reactor operation. During normal reactor operations, a power increase can be brought about with the help of methods similar to those employed during the BPGCR startup, i.e., a combination of external positive reactivity addition (e.g., by means of control rods or drums), plus a change in the core inlet pressure. This method is discussed in detail in the results section of Chapter 3. To see whether the delayed neutron circulation effect is influencing the system power behavior significantly during reactivity transients, two additional transients are studied; in the first one, the external loop residence time of the fuel is increased from 0.1 to 0.3 seconds, and the same positive external reactivity step of magnitude 0.0029 Ak/k is applied to the system operating at a steady state power of 100 MWt. In the second one, the external loop circulation time is increased from 0.1 seconds to 0.5 seconds and an external positive reactivity step of magnitude 0.0029 Ak/k is applied to the system (the core fuel residence time is kept fixed at 0.08 seconds in both cases). A change in the external loop circulation time affects the BPGCR
PAGE 122
114 dynamics by changing the fraction of delayed neutron precursors that are reintroduced into the core (without deca 5 dng) after traversing the external loop. Figures 4.5 and 4.6, respectively, show the system power behavior for the above two cases. For the purpose of comparison, the results from the transients shown in Figs. 4.1, 4.5, and 4.6 are combined and are shown together in Fig. 4.7. From Fig. 4.7, it can be seen that the variations in the external loop circulation time have very little effect on the initial oscillatory behavior of the system (the difference in amplitudes of the maximum power attained by the system, immediately following the external step reactivity insertion, is less than 5% for the three cases shown) but, the changes in the external loop circulation time do significantly affect the final power attained by the system at the end of the transient (there is a difference of roughly 25% in the final power attained by the BPGCR among the three cases shown in Fig. 4.7). The initial oscillatory behavior is primarily dictated by the strong, prompt fuel mass feedback effect, while the final power is dictated by the much slower delayed neutron feedback effect due to the fuel circulation. The final powers are higher for the cases with longer external loop circulation times (Ti = 0,5 second and Ti = 0.3 second) compared to the case with the shorter external loop circulation time (li = 0.1 second) because of the fact that, although the same (magnitude) external reactivity step in absolute units is applied to the system in all cases, for the longer
PAGE 123
1.3E+8 115 raajsiis rs^ox Figure 45. BPGCR Thermal Power Vs. Time After an External Reactivity Insertion Effect of Varying External Loop Circulation Time Case 1
PAGE 124
Initial Core Power = 100 MWt 116 u 9i m fi O t: V (A a HH M Â•5: Â•M U ed V 0^ Oi V M (0 73 e V H 4) g H Figure 46. BPGCR Thermal Power Vs. Time After an External Reactivity Insertion Effect of Varying External Lx)op Circulation Time Case II
PAGE 125
1.2E+84 Initial Core Power = 100 MWt 117 Figure 47. BPGCR Thermal Power Vs. Time After an External Reactivity Insertion Effect of Varying External Loop Circulation Time
PAGE 126
118 external loop circulation time cases (i.e., with Ti = 0.5 second and li = 0.3 second), the external reactivity step corresponds to more than a dollar in relative units. For a circulating fuel reactor, the magnitude of relative reactivity that corresponds to a dollar is given by ' a^', where 'a' is defined as the delayed neutron attenuation factor. The factor 'a' depends on the core residence time as well as on the external loop circulation time of the fuel (see Eq. 4.4 for an expression for 'a'). The magnitude of the applied external reactivity step corresponds to exactly a dollar for the case with Xi = 0.1 second (and Tc = 0.08 second). This difference in reactivity insertion, measured in relative units, is why the amplitude of oscillations for the case with Xj = 0.3 second is slightly higher compared to the case with Xi = 0. 1 second, and the amplitude of oscillations for the case with Xi = 0.5 second is slightly higher compared to the cases with Xi = 0. 1 second and Xi = 0.3 second. To help determine the error associated with performing a dynamic analysis of the BPGCR with only one delayed neutron group (as opposed to performing the analysis with all six delayed neutron groups), a twenty cent (20
PAGE 127
1.2E+8 119 rB?ox One Delayed Neutron Group Case
PAGE 128
1.2E+8 120 o d o 00 o d o q d o d (s;;bm) JSMOd ma^sXs ib:jox o V Â•21 a o t: V (A a o ci Â«i (H a M (0 cd g M X u a B H Figure 49. BPGCR Thermal Power Vs. Time After an External Reactivity Insertion Six Delayed Neutron Group Case
PAGE 129
121 that, for the examined system, neither the power behavior immediately folloAving the perturbation, nor the final steady state power attained by the system are significantly affected by inclusion of all six delayed neutron groups. Comparing Figs. 4.8 (one delayed neutron group case for 20t positive reactivity step) and 4.3 (one delayed neutron group case for $1 negative reactivity step), it can be seen that the final steady state power attained by the BPGCR is almost the same. This behaviorÂ— as mentioned beforeÂ— is very different from that of a conventional, solid core, stationary fuel reactor. In a conventional reactor, the fuel mass within the reactor remains constant and, if a positive (or negative) external reactivity step is applied to such a system operating at steady state, the net system reactivity will remain positive (or negative), and the system power will correspondingly increase (or decrease) after the step reactivity insertion. For the BPGCR, the fuel mass within the core does not remain constant during the transient. Immediately following the reactivity perturbation, the amount of fuel mass within the core undergoes oscillations before quickly coming to an equilibrium. The magnitude of the equilibrium internal reactivity (associated with the fuel mass within the core) is such that the total reactivity (internal plus the externally applied) will bring the system to a new steady state power level. Note that because of the loss of delayed neutron precursors from fuel circulation, the magnitude of the reactivity needed for steady state power is greater than zero.
PAGE 130
122 Since the fuel mass feedback effect is dictating the system power behavior during the transients described above, it is decided to see the effect of changing the magnitude of the fuel mass feedback coefficient of reactivity on the system behavior. Under normal conditions, the value of this coefficient changes with time depending on the operating conditions (i.e., the fuel mass within the core dictates the magnitude of this coefficient). To see the effect of changing this coefficient, within DYNAM, the magnitude of this coefficient (Ap/AMf) is artificially changed, and is fixed at two separate values of 0.001 and 0.0005 (Ak/k per Kg of fuel), respectively (during normal operating conditions, the magnitude of this coefficient is of the order of 0,035). The results obtained are given in Figs. 4.10 and 4.11, respectively. As the value of the fuel mass feddback coefficient is reduced, the oscillations previously exhibited by the system (Figs. 4.1 through 4.6), disappear completely. This is an expected result since, we are essentially decreasing the negative feedback mechanism available to the system. The fuel mass is no longer as strong a feedback mechanism as it was before, and the variations in core pressure and core temperature do not affect the system reactivity as much as they did previously. But if we reduce the fuel mass feedback coefficient too much, the system becomes unstable and the power starts increasing without limit after reactivity perturbations (see Fig. 4,11).
PAGE 131
3.0E+8 123 jsmod ma^s^s I^^ox Figure 410. BPGCR Termal Power Vs. Time After an External Reactivity Insertion Effect of Varying the Fuel Mass Feedback Coefficient Case I
PAGE 132
6.0E+8 124 nZBl^SiCs I^^OX Figure 411, BPGCR Thermal Power Vs. Time After an External Reactivity Insertion Effect of Varying the Fuel Mass Feedback Coefficient Case II
PAGE 133
125 At the other extreme, too large a value for this coefficient will also make the system unstable. For example, results for the case where the value of this coefficient is increased from its typical value of 0.035 to a value of 0.1, are shown in Fig. 4.12. The oscillations now take more time to die down compared to previous cases, and it is found that the BPGCR will break into undamped oscillations if the value of this coefficient is increased beyond 0.2. An optimum value of this coefficient can be selected with the help of figures like the ones shown and the system design can be modified so as to operate in a regime of the fuel mass versus reactivity curve so as to obtain the desired value of afmThe coefficient afm is obtained from the slope of the keff versus core fuel mass curve. The value of ttfm can, thus, be changed by operating at a different fuel mass loading. In order to maintain the needed keff, this will require some form of reactivity compensation. This can be achieved with the help of external control rods or drums located in the external moderator/reflector region or by adjusting the gas fuel enrichment or by adding a burnable poison to the fuel gas mixture. Finally, to the see the effect of changing the inlet pressure to the core, at the same time that the one dollar reactivity step insertion is applied to the core, the value of the core inlet stagnation pressure is
PAGE 134
1.2E+8 126 (s;;bm) J3MOJ uza:^SiCs I^^ox Figure 412. BPGCR Thermal Power Vs. Time After an External Reactivity Insertion Effect of Varying the Fuel Mass Feedback Coefficient Case III
PAGE 135
127 raised to a much higher value of roughly 38 atm. from its typical value of about 33 atm. This core inlet stagnation pressure is kept at a constant value for this as well as for all transients described previously. Since the fuel mass feedback mechanism is dictating the system power behavior during the transient, we expect any effects that affect the fuel mass inside the core to affect the power behavior also. This can be seen from Fig. 4.13, where the total core power is plotted as a function of time after the initiation of the transient. Since, now we have a higher core inlet pressure, more fuel is forced into the core and, hence, we get higher power oscillations (due to higher internal reactivity oscillations). Also, the frequency of oscillations have increased now compared to the previous cases. This is an expected result since the frequency of these oscillations depends on the fuel transit time through the core which depends on the fuel inflow and outflow rates which are in turn affected by the higher core inlet stagnation pressure. The relative magnitude of the core inlet pressure with respect to the average core pressure affects the core fuel loading and, hence, the final steady state power level attained by the system after an external reactivity perturbation. The magnitude and sign of the externally applied reactivity influences only the initial power oscillations that take place immediately following this perturbation. The effect of varying the core inlet stagnation pressure can be seen from Figs. 4.14 through 4.16, where the BPGCR power behavior
PAGE 136
Initial Core Power =100 MWt $ 1 External Positive Reactivity Step 128 jaMOj ni3:(s^s I^^ox o V (A Qi 9i (/i Â•a e V Â•M X a IH V V B Figure 413. BPGCR Thermal Power Vs. Time After an External Reactivity Insertion Effect of Increasing the Core Inlet Stagnation Pressure
PAGE 137
Core Average Pressure Before Reactivity Insertion = 33 atm Constant Core Inlet Pressure = 34.1 atm 129 o in 00 + O jaMOj I^^ox Figure 414. BPGCR Thermal Power Vs. Time After an External Reactivity Insertion Effect of Varying the Core Inlet Stagnation Pressure Case I
PAGE 138
1.2E+8 130 o j3MO({ xna^s^s I^^ox u V m a Â•M (/) Â•a e 0> M M U Â»4 V Figure 415. BPGCR Thermal Power Vs. Time After an External Reactivity Insertion Effect of Varying the Core Inlet Stagnation Pressure Case II
PAGE 139
Core Pressure 131 o (s^;bm) J3AOd uzai^s^s I^^ox Figure 416. BPGCR Thermal Power Vs. Time After an External Reactivity Insertion Effect of Varying the Core Inlet Stagnation Pressure Case III
PAGE 140
132 following a one dollar positive reactivity step is shown for three different values of constant core inlet stagnation pressure. It can be seen from the figures that, as the value of the core inlet stagnation pressure is increased (relative to the average core pressure), the final core power, as well as the amplitude of the initial power oscillations, also increase. For all the reactivity perturbations examined so far, the core inlet stagnation pressure is kept constant during the transient. It is established that the BPGCR is an extremely stable configuration under this condition. By keeping the inlet pressure constant during the transient, we are, in effect, providing an external controlling mechanism. To see whether the BPGCR is stable even without this external action (i.e., to see whether the system is inherently stable or not), a one dollar positive external reactivity step is applied to the system and the core inlet pressure is allowed to vary. The core inlet pressure at any given time, t, is set equal to the average core pressure which had existed at a time Xi back, i.e., it is set to the average core pressure at time (tTi), where Xi is the external loop circulation time. The resulting BPGCR power behavior in Fig. 4.17, shows that the BPGCR thermal power drops more than two orders of magnitude in a matter of few seconds (note that the Yaxis of Fig. 4.17 has a log scale). So, it can be seen that when the core inlet pressure is allowed to float, the BPGCR becomes even more stable against accidental reactivity
PAGE 141
6+aOT 133 jamoj ma^s^s I^^ox Figure 417, BPGCR Thermal Power Vs. Time After an External Reactivity Insertion Variable Core Inlet Stagnation Pressure Case
PAGE 142
134 insertions. For normal reactor operation, this power behavior is much less desirable than for the constant core inlet pressure case because of the large drop in power Avithin a short period of time, as well as the fact that the power level continues to drop without coming to a new steady state power level (not shown in Fig. 4.17). The above described method of analyzing the transient response (stability) of a reactor by solving simultaneously the set of nonlinear system equations can rapidly become impractical since the complete reactor system is quite complex, and usually many equations are required to accurately describe the behavior of the system, and many transients for many different values of system parameters need to be examined to map the entire operating range. Solving these nonlinear system equations can be both time consuming as well as tedious. A second, and less complicated, approach to studying the stability is to linearize the feedback terms in the system equations and use the welldeveloped methods of linear feedback theory for stability analysis. Some important conclusions about a nonlinear system can usually be made by stud 5 ang its associated linear approximation. In particular, the stability of a nonlinear system for small perturbations may often be deduced by examining the stability of an associated linear system. This type of analysis is particularly useful for obtaining a good understanding of the general nature of the reactor behavior. Effects of changes in input reactor parameters can be readily identified, and
PAGE 143
135 regions of stability can be mapped for various ranges of reactor parameters. Such a linear stability analysis of the BPGCR is undertaken in the following section. The methods of analysis include the familiar frequency response methods (Bode plots, Nyquist plots, etc.) and pole configuration (root locus) methods. Details of these methods are not given in this dissertation, but can be found in a variety of system dynamics textbooks including Ref, 26. BPGCR Linear Stability Analysis Dynamically the most important characteristics of the gas core reactor being investigated in this study are that the fuel circulates and that the fuel is highly compressible. The fuel circulation acts to reduce the effective delayed neutron fraction (in general, the lower the value of delayed neutron fraction, the more difficult it becomes to control the reactorÂ— a reason why plutonium fueled reactors are more difficult to control than uranium fueled ones). From a control standpoint, the delayed neutron fuel circulation effect acts as a positive feedback (i.e., if the reactor power goes up, the delayed neutrons reintroduced into the core after going through the external loop will also go up, and this tends to drive the reactor power up further) while the gas density effect is a negative feedback. These are two opposite (in sign) but powerful feedback control mechanisms. Fuel circulation also acts to reduce the rate of fuel temperature change during a power transient. For a linearized system, the transfer function is defined as the Laplace transform of a selected output variable divided by the Laplace
PAGE 144
136 transform of a selected input variable. Its value is a complex number describing both amplitude and phase relationships. For a reactivity input, the reactor transfer function is the ratio of the Laplace transform of the steady state flux response to the Laplace transform of the reactivity input, i.e., for an input reactivity p(t) = po Sin (Ot, the flux will, eventually (after the transients have died down), be given by (t) = 00 sin (cot + 0) and the transfer function is given by G Oco) = ^ exp (j0) PO where '0o' is the amplitude of the of the flux response, 'po' is the amplitude of the reactivity input, and '0' is the phase shift. The open loop transfer function is the response of the reactor in the absence of any feedback (temperature coefficients of reactivity, etc.). This transfer function is characteristic of a reactor operating at zero power where the flux level is so low that changes in the flux level do not significantly affect the temperature or other feedback variables. The closed loop transfer function is the response of the reactor at power, including all feedback mechanisms, whereby the power level affects the reactivity of the reactor. References 26 and 29 contain more details about transfer functions, and how they can be derived from the system equations, and how they can be employed to obtain information about the system stability.
PAGE 145
137 If we consider an elementary critical reactor (i.e., one with no feedback and no fuel circulation), the transfer function of this reactor can be easily obtained from the neutron kinetics equations (see Ref. 26), and is given below for the case of one delayed neutron group. G(s) = ^ (s + A.) A s (s+ p/A) (41) where no = equilibrium value of the neutron density (#/cc), s = Laplace transform variable, X = average value of decay constant for the delayed neutron precursors (1/sec), A = neutron generation time (sec), and 3 = total delayed neutron fraction. If we analyze this transfer function of an elementary reactor using any of the standard linear stability analysis tools, we will find that this reactor is unstable. For example. Fig. 4.18 shows the Bode plot of this transfer function. A Bode plot is a simple plot of the magnitude and phase angle of the transfer function versus frequency (often. Bode plots are used to display experimental results in a form from which an approximate transfer function can readily be derived). Figure 4.18 (in which only the amplitude plot is given for the sake of illustration) shows that the reactor becomes unstable at zero frequency, i.e,, the reactor has infinite gain at zero frequency. Physically, this means that
PAGE 146
l.OOE+13 138 Figure 418. Bode Plot for Reactivity Transfer Function of an Elementary, Zero Power, NonCirculating Fuel Reactor
PAGE 147
139 when a small amount of positive reactivity is inserted into a critical reactor, its power level rises indefinitely. For this Bode plot, 'no' is taken to be equal to 1.0 E+8, and the rest of the parameters have the same value as were used in DYNAM (i.e., a delayed neutron fraction of 0.0065, a neutron generation time of 6.0E4 sec., and an average precursor decay constant of 0.08 seci). In a real reactor, whenever operating at power, there are always feedback effects present which tend to alter the response of the reactor to external inputs or perturbations. When feedback effects are present, the elementary reactor transfer function is modified as shown in Fig. 4.19, and the overall reactor transfer function is given by Gtot(s) = G(s) / (1 + G(s)H(s)) (42) where G(s) = elementary zero power reactor transfer function, and H(s) = feedback transfer function If we now include the effect of fuel circulation on the elementary zero power reactor model, we have to modify the elementary reactor transfer function given above. The circulating delayed neutron precursors tend to act like a delayed positive feedback on the reactor transfer function. To understand this phenomenon, consider a circulating fuel reactor operating at steady state. Of all the delayed neutrons precursors produced in the core, only a certain fraction will decay within the core before leaving the core. Some of the precursors
PAGE 148
140 u rt T3 (U (U (Zh c o 3 u c 3 b u CO tl H ii o M u (U Ci O) d) Ih 3 130
PAGE 149
141 will decay while they are still in the external loop and, hence, will not contribute to the core neutron population. The ones that survive without decaying while traversing the external loop will reenter the core and then can decay. When the reactor is operating at steady state, a dynamic equilibrium is established for the delayed neutron population within the core. Now, if for some reason, the core power goes up, the delayed neutron precursor density will also increase. This increase in the precursor density manifests itself as an increase in the fraction of precursors that are reintroduced into the core after some delay (depending on the external loop circulation time) and, hence, will tend to drive the power further up. The reverse of this process occurs when the power goes down. Using the equations given in Ch.3 (Eq. 3.12) for the precursor density in a circulating fuel reactor, one can derive the transfer function for a zero power, circulating fuel reactor (see Ref. 29 for details of the derivation). It is given by the following expression: G(s) = MÂ£) = s2I + s ja + + apj + XP(al) + where C 2 = 1 e^'^i e Â®^i Xc = mean core fuel residence time (sec).
PAGE 150
142 Xi = mean external loop eirculation time (sec), X = average decay constant for delayed neutron precursors (seci), = total delayed neutron fraction, and a = delayed neutron attenuation factor (see discussion below) In the above equation, the quantity "a" is called the delayed neutron attenuation factor. Its effect is to reduce the delayed neutron fraction in a circulating fuel reactor by a factor "a". We can intuitively infer that the effective delayed neutron fraction in a circulating fuel reactor has to be less than that for a noncirculating fuel case. Also, intuitively, we can say that this factor "a" has to depend on the relative times the precursors spend within the core and in the external loop. For example, if the external loop circulation time is large, a greater fraction of delayed neutron precursors will decay while in the external loop. We can formally derive the expression for "a", to obtain a = ^ Xxc + ie^^i (44) (See Ref. 30. for details of this derivation.) For the gas core reactor being investigated, in the program DYNAM, the external loop circulation time is generally fixed at 0. 1 second. The core residence time of the fuel depends on the velocity of fuel within the core but is generally in the range of a few tens of milliseconds (2040) corresponding to fuel velocities of the order of 50 to 100 m/sec (assuming a core length of about 2 meters). The Bode
PAGE 151
143 plot (magnitude) for this zero power, circulating fuel reactor transfer function is given in Fig, 4.20 (a 5th order Fade approximation is used to represent the delay terms of Equation (4,3)). Comparing this figure with Fig. 4. 18 (which shows the transfer function for a noncirculating fuel reactor), it can be seen that the gain is higher (i.e., the system is less stable) at lower frequencies for the circulating fuel reactorÂ— an expected result due to the lower effective value of the delayed neutron fraction for the circulating fuel reactor. Inherent Stability of the BPGCR Attached to an External Loop In the previous section, the reactivity transfer function of the BPGCR operating at zero power was given. However, in a reactor operating at power, feedback effects are always present, and an external loop is always attached to the reactor. The ultimate concern is the stability of the complete system, i.e., we must include the effects from within the reactor, as well as the effects of the balanceofplant (B.O.P), to check for the overall stability. These internal feedbacks and the B.O.P can significantly alter the response of a zero power reactor model given in the previous section. From an elementary transfer function standpoint, the complete system can be considered as shown in Fig. 4.21. Here all the possible feedback effects that might affect the stability of the system are represented as two distinct loops. Each one of these two loops can be further broken down to include component transfer functions.
PAGE 152
CT+aOOT 144 dpn^lTZ^BIfQ Figure 420. Bode Plot for Reactivity Transfer Function of a Zero Power, Circulating Fuel Reactor
PAGE 153
145 Figure 42 1 . Simplified Schematic Diagram of the Complete BPGCR System
PAGE 154
146 The internal loop includes the various power coefficients of reactivity, i.e., functions which describe how the reactor power influences the reactivity effects due to changes in the fuel temperature, the moderator temperature, the fuel mass, etc.. The external loop represents the effect of components like the heat exchanger, the radiator, the compressor, etc., on the performance of the BPGCR. These effects include the heat transfer from the fuel/ coolant to the heat exchanger (attenuation of power), the time involved in the process of circulating the fuel through the external loop (time delay), etc. The combination of these external loop transfer functions, in conjunction with the zero power reactor transfer function and its associated internal feedback functions, determine the stability of the overall system. External Loop As mentioned before, the external loop is associated with the B.O.P. For the BPGCR, the fuel and the working fluid (coolant) are one and the same. The gaseous fuel/ working fluid mixture leaving the reactor will remove most of the heat generated within the core. Under steady state conditions, this heat is given up by the fuel gas mixture while going through the external loop and the fuel reenters the core at a lower temperature after some delay. The amount of attenuation (in temperature) and delay the fuel gas mixture suffers while going through the external loop depends on the components of the external loop. To see how the external loop can affect the performance of the
PAGE 155
147 reactor, consider a sudden increase in the power generation within the core. This results in an increase in the temperature of the coolant (fuel) exiting from the core and, hence, that of the fuel gas mixture reentering the core. An increase in these two temperatures will affect the average fuel temperature and fuel gas mixture density within the core which, in turn, can influence all the internal feedback loops of the reactor. To derive the transfer functions for the external loop, first the components of the loop have to be selected and then the equations describing their performance have to be modelled. Since the BPGCR is only in a conceptual stage of development, there is a lot of flexibility in selecting the various components of the B.O.P. and their operating conditions (for example, the nature of the working fluid in the secondary side of the radiator/heat exchanger or the average secondary fluid temperature, etc.). For the purpose of this analysis, a simple external loop as shown in Fig. 4.22 is selected. It consists of a turbine (in the final design there may not be a turbine as such but some sort of an energy conversion device such as an MHD generator will be thereÂ— which is represented here as a turbine), a heat exchanger /radiator combination for rejecting the waste heat, and finally a compressor for pumping the fuel back to the core at the appropriate pressure. The equations describing each of the above components are given in Appendix B. Since this is a linear stability analysis, only small
PAGE 156
148
PAGE 157
149 variations around some mean operating conditions are considered, i.e., all parameters are treated as perturbation quantities or deviations from their steady state values. The transfer function representation of the external loop can be represented as shown in Fig. 4.23 by using the equations given in Appendix B. As shown in the figure, the total time delay in the external loop is handled by two lumped delays (G 6 )Â— one for the transfer of fuel gas mixture from the core to the heat exchanger and the other for the fuel gas mixture transfer from the heat exchanger back to the core. In Fig. 4.23, the transfer function G5, represents the combined effect of the turbine, the heat exchanger/radiator, and the compressor. Gl, G2, G3 and G4 represent the transfer functions that describe how the core outlet temperature (T 3 ') is affected by the core inlet temperature (T 2 ), core thermal power (N). core inlet fuel mass flow rate (m 2 ), and the core exit fuel mass flow rate (m 3 ), respectively. Appendix B gives the contents of the transfer function boxes shown in Fig. 4.23. Internal Feedback Loon The structure of this loop for the BPGCR is significantly different form that of a conventional, solid core reactor. Again, the main reason for this difference is the gaseous nature of the fuel plus the fact that BPGCR is an open flow system. For the BPGCRÂ— unlike for a conventional solid core PWR, the most significant component of the reactor internal feedback loop is the fuel mass (or density) coefficient
PAGE 158
150 Figure. 423. Transfer Function Representation of the BPGCR External Loop
PAGE 159
151 of reactivity. In a conventional solid core reactor, the fuel mass within the reactor remains constant whereas in the BPGCR, the gaseous fuel mass is a strong function of the operating conditions (core average temperature, pressure, etc.) and, hence, constitutes a strong feedback mechanism. The other components of the internal feedback loop like the fuel temperature (Doppler) coefficient of reactivity and the moderator temperature coefficient of reactivity are of much less importance compared to the fuel mass coefficient of reactivity. The equations describing the internal feedback loop are given in Appendix B. As a first step in the analysis of the internal feedback loop, it can be represented as shown in Fig. 4.24. This figure shows the three major feedback loops that can influence the stability of the BPGCR. The core inlet and outlet temperatures combine to provide a core average fuel gas mixture temperature. Changes in this average fuel gas mixture temperature directly gives rise to one component of the internal feedback loop, viz., the fuel temperature (Doppler) coefficient of reactivity. The fuel gas mixture temperature also indirectly influences the moderator temperature and the fuel mass (density) components of the internal reactivity feedback loop as shown in the figure. Formulations for the fuel temperature feedback and the moderator temperature feedback parts of the internal feedback loop are relatively straightforward. The average core fuel temperature is obtained as a simple average of the core inlet and outlet temperatures.
PAGE 160
G7 152 Figure 424. Transfer Function Representation of the BPGCR Internal Reactivity Feedbacks
PAGE 161
153 and the fuel temperature coefficient of reactivity is estimated using 123group XSDRNPM (Sn transport theory) calculations (see Chapter 3. Tables 3.5 and 3.6). For the moderator temperature coefficient component, the equations relating the variation of the average moderator temperature as a function of the variation of the average fuel temperature are given in Appendix B. The moderator temperature coefficient of reactivity is also estimated using 123group XSDRNPM calculations (see Chapter 3, Table 3.7). The analysis of the gaseous fuel mass (or density) component of the internal feedback loop is more involved and has to be performed in steps. The fuel mass inside the core is affected by a number of variables including the core pressure, core temperature, core inlet stagnation pressure, fuel mass inflow and outflow, etc.. The complete set of equations interrelating all these variables are given in Appendix B. Rearrangement of these equations and some algebraic manipulations (see Appendix B), yield the block diagram given in Fig. 4.25which shows the relationship among the involved variables. For the initial stability analysis, we will consider only the fuel mass component of the internal feedback loop. This is the strongest of the three internal feedback paths. The fuel mass feedback effect can be added to Fig. 4.25 since it gives the variation of fuel mass within the BPGCR. Adding the fuel mass feedback coefficient to this fuel mass variation, yields the variation of that component of internal reactivity which is due to the fuel mass variation. Adding a zero power
PAGE 162
G7 154 CO < + M H CO Figure 425. Transfer Function Representation of the BPGCR Fuel Mass Feedback
PAGE 163
155 circulating fuel reactivity transfer function to this internal reactivity variation gives the thermal power variation within the BPGCR. This is shown in Fig. 4.26Â— where, in addition to the fuel mass feedback component, the transfer functions which represent the influence of core inlet temperature, core thermal power, core inlet and exit fuel mass flow rates on the core exit temperature (see Fig. 4.23), are also included. Figure 4.26 now represents the complete BPGCR except for the balanceofplant. Before the stability analysis of the system can be performed, the mass feedback branch of the internal feedback loop given above has to be further reduced. After reducing the block diagram given in Fig. 4.26, we can combine it with the balanceofplant loop given earlier (Fig. 4.23) to obtain a simplifled representation of the complete BPGCR system as shown in Fig. 4.27. This is the model that is used for the stability analysis of the system. The linear stability analysis of the BPGCR is performed using a computer program called CC. It is a personal computer based commercial software package which is capable of performing both classical control theory analysis, as well as state variable analysis. For performing a classical stability analysis (frequency response, root locus, etc.) using CC, the complete system given in Fig. 4.27 has to be reduced to a single transfer function, i.e., a single block representing the relationship between an input and an output as shown below:
PAGE 164
G7 156 CO Figure 426. Transfer Function Representation of the Complete BPGCR with Fuel Mass Feedback Excluding the BalanceofPlant
PAGE 165
G21 157 CO Figure 427. Simplified Transfer Function Representation of the Complete BPGCR
PAGE 166
158 Input Output To get the system block diagram into such a form, first the external loop can be eliminated. To do this, we have to derive two transfer functions: one representing the variation of the core outlet temperature, T 3 ', with respect to the variations in the core power ( 3 T 3 / 3N), and a second one representing the variations of the core inlet temperature, T 2 ', with respect to the variations in the core power OT 2 /3N). These two transfer functions can be derived easily from the external loop representation given in Fig. 4.27 (the derivation is given in Appendix B). With the help of these two transfer functions, the block diagram of Fig. 4.27 can be reduced to the one shown in Fig. 4.28. Fig. 4.29 shows the same overall transfer function block diagram in a more conventional manner, where an external reactivity input is also shown. This external reactivity could, for example, represent the control rods or an adjustable external reflector, or some type of an accidental reactivity insertion. The total system reactivity is given by the sum of this external reactivity and the internal reactivity. The block shown in Fig. 4.29 can now be used for obtaining the overall system transfer function, 9N / 3pext. as
PAGE 167
159 Figure 428. Reduced Transfer Function Representation of the Complete BPGCR Figure 429. Reduced Transfer Function Representation of the Complete BPGCR in Standard Form
PAGE 168
dN ^ ZPR dpext 1ZPR am (G22 G24 + G21 G23) 160 (45) where ZPR is the zero power, circulating fuel reactor transfer function. Linear Stability Analysis of the Overall System Transfer Function Equation 4,5 is of the familiar closed loop transfer function form Y(s) = G(s) / (1 + G(s)H(s)) where G(s) represents the feed forward transfer function and H(s) represents the feed back transfer function. The transient response of the system is characterized by the poles of Y(s)Â— which are the roots of the characteristic equation 1 + G(s)H(s) = 0 (46) or, in our case, 1 + (ZPR ttM (G22 G24 + G21 G23)) = 0. (47) Linear stability is assured if all the poles of the system transfer function are in the left half plane, i.e., if roots of the characteristic equation are negative or have negative real parts. Using the program CC, the roots of the characteristic equation, Eq. 4,7, are obtained and are given below: Zeros of 1 ZPR Qm (G22 G24 + G21 G23) = 0 2.1850 E4 7.8550 E2
PAGE 169
161 1.830 + 0.5463 Â± j 4.0518 50.4184 5.3712 Â±j 47.4075 3.6892 Â± j 54.2463 338.486 26.888 Â±J 13 1.7 18 and, 900.0180 The above values for the roots of the characteristic equation are obtained for the following BPGCR operating conditions: Core Thermal Power = 100 MWt, Core Average Fuel Temperature = 1500 K, Average Beryllium Temperature = 800 K. Average Radiator Temperature = 1200 K, Core Inlet Stagnation Pressure = 31 atm. Fuel Mass Flow Rate = 140 Kg/s, Fuel Mass Feedback Coefficient = 0.035 Ak/k per Kg, Fuel External Loop Circulation Time = 0.1 sec, and Fuel Core Residence Time = 0.08 sec It can be immediately seen that, because of the positive root (complex conjugate pair) on the right half plane, the system is unstable, i.e,, as time progresses after a reactivity perturbation, the disturbance (system output) grows without limit.
PAGE 170
162 Earlier in this chapter, when the nonlinear analysis of the system was performed using DYNAM, it was found that the fuel mass feedback is a strong feedback mechanism, and the system can be rendered stable or unstable depending on the value of this feedback coefficient. It is decided to see whether the fuel mass feedback has the same effect in a linearized model also. To verify this, the value of the fuel mass coefficient is reduced from its typical value of 0.035 to a value of 0.007 (a factor of two reduction). The roots of the characteristic equation are again determined using CC and are given below; 2.1850 E4 7.2280 E2 1.830 0.4687 Â±j 3.336 52.177 4.4060 Â±j 26.2160 3.6810 Â±j 54.2480 338.486 26.888 Â± j 131.717 and, 898.1870 It can be seen that the system is now stable for this value of the fuel mass feedback coefficient. It is sometimes much easier to visualize the system stability by looking at the system response in time domain. Using CC, the time
PAGE 171
163 response of the BPGCR (linearized model) is studied for various values of the fuel mass feedback coefficient, and are given in Figs. 4.30 through 4.33. Since these results are obtained using a linearized model only a small reactivity step (10 cents) is applied (large reactivity perturbations will render the linearized model invalid). Comparing these figures with Figs. 4.10 through 4.12, it can be seen that the linearized model predicts the same type of BPGCR response as the nonlinear model, i.e., the system becomes unstable for both very large and very small values of the fuel mass feedback coefficient. There is a range of values of this coefficient where the system is stable to external reactivity perturbations. The actual values of the fuel mass feedback coefficient (both high and low) at which the system becomes unstable as predicted by the two models (linear and nonlinear), however, are different. Now that it is determined that the system is stableÂ— at least for certain regimes of the fuel mass feedback coefficient, the next step is to determine whether this is true for all values of the system gain (reactor power). This can be determined from a root locus plot of the transfer function of the system. The root locus method is based on the relationship between the poles and zeros of the closed loop transfer function and those of the open loop transfer function. The quantity Â’G(s)H(s)' in Equation 4.6 is called the open loop transfer function, and it plays a very important part in the stability analysis of the system (see Ref. 26 for details of the root locus method). Basically, in the root
PAGE 172
1.50E+7 164 Figure 430. Linear Analysis Prediction of the BPGCR Thermal Power Vs. Time After an External Reactivity Insertion Case I
PAGE 173
2.0E+7 165 Figure 431. Linear Analysis Prediction of the BPGCR Thermal Power Vs. Time After an External Reactivity Insertion Case II
PAGE 174
2.0E+7 166 (s^;bm) J3AiO(i uza^si^s ps^ox Figure 432. Linear Analysis Prediction of the BPGCR Thermal Power Vs. Time After an External Reactivity Insertion Case III
PAGE 175
4.0E+7 167 j3Avod I^^ox Figure 433. Linear Analysis Prediction of the BPGCR Thermal Power Vs. Time after an External Reactivity Insertion Case IV
PAGE 176
168 locus method, the behavior of the roots of the characteristic equation (the loci) are investigated as a function of the system gain. A necessary and sufficient condition for linear stability is that the roots all have negative real parts. If we write the characteristic equation in the form 1 + K G(s) H(s) = 0. then the root locus method gives the loci of its zeros as K varies from zero to infinity (the gain K include the reactor power and the fuel mass feedback coefficient as a part of it). The open loop transfer function can be expressed as a ratio of two factored polynomials rfci Rfci K (sri) (sr 2 ) . . . (srj (spi) (sp 2 ) . . . (spm) (48) where 'ri' are called the zeros and 'pi' are called the poles. Fig. 4.34 shows the root locus plot of the overall BPGCR transfer function for the same operating conditions given earlier. The 'X's in the figure denotes the poles and the 'O's denotes the zeros of the open loop transfer function. There are three asymptotes, approaching infinity along 180, +60, and 60 degrees (Ref. 26 gives the rules for constructing a root locus plot). In the figure, for the sake of clarity, the real and imaginary axes are not drawn to scale. It is meant for the purpose of illustrating the behavior of the roots of the characteristic equation only. From the figure we can see that four loci cross over to the right hand side (positive X coordinate side) as the system gain increases
PAGE 177
Figure Not to Scale 169 Figure 434. Root Locus Plot of the BPGCR with Fuel Mass Feedback
PAGE 178
170 from zero to infinity. The system gain at the cross over points are also shown in Fig. 4.34. It can be seen from the figure that the limiting value of the system gain is about 2.5, i.e., the system is unstable for all values of gain, K above 2.5. Also, from the root locus plot, as the system gain increases, we can expect the oscillatory behavior of the system output to increase. So, a linear stability analysis also indicates that the BPGCR system analyzed in this dissertation is a stable configuration for a range of operating conditions of the fuel mass feedback coefficient.
PAGE 179
CHAPTER 5 CONCLUSIONS, RECOMMENDATIONS AND AREAS OF FURTHER RESEARCH Introduction The d 3 niamic analysis of a conceptual gas core reactor performed in this dissertation indicates that the Burst Power Gas Core Reactor system or, the BPGCR system, is capable of rapid startup and very stable operation at high power levels in a space environment. In this dissertation, the BPGCR is primarily considered for operation in a space environment, even though it is equally capable of terrestrial operation. It is found that the combination of the gaseous nature of the fuel, plus the circulating fuel configuration introduces unique feedback effects which tend to alter significantly the dynamic response of the BPGCR compared to that of conventional, solid core systems. These unique feedback effects are also found to make the system very stable during operation at high power levels. The BPGCR considered in this dissertation is designed as an "open flow" system with nozzles, rather than valves or orifices, at the core inlet and at the core exit. In addition to the simplicity of design, such a configuration can offer some added benefits from a materials point of view when the BPGCR is 171
PAGE 180
172 operated using highly corrosive gaseous fuel forms like UFe or UF 4 . The analysis performed in this dissertation indicates that this BPGCR design, with only nozzles at the inlet and exit, is an inherently stable configuration during high power mode of operation. Results Transition Analysis Since the BPGCR is a novel system, no previous analysis exists on the dynamics of such a system, nor is any information available on the system behavior during rapid startups. In this dissertation, a computer program (DYNAM) is developed to perform this analysis. The computer program DYNAM solves simultaneously the coupled neutronics, heat transfer, and thermodynamic equations of the BPGCR attached to a balanceofplant. Since no similar program exists, various sections of DYNAM are benchmarked against standard computer codes, or, in some cases, against analytical results. DYNAM can be a very valuable tool for analyzing the startup of the BPGCR and also for further refining the system design. Using DYNAM. it is established that the BPGCR is capable of going safely and in a predictable manner from a low power mode of operation to a high power mode (a few 100s of MWt) of operation in a short period of time (of the order of 100 seconds). Additionally, it is found that only simple external controlling action is needed to control the system during this rapid transition. In this dissertation, the
PAGE 181
173 method considered is to control the system by controlling the core inlet pressure. Different startup scenarios (different types of reactivity insertions) are considered in this dissertation, DYNAM can predict the behavior of important system parameters (like core pressure, temperature, fuel mass flow rates, etc.) during any of these startup scenarios. An advance knowledge of the behavior of these parameters can then enable one to select desirable transition paths, i.e., paths along which the crucial parameters of the BPGCR are within satisfactory bounds. DYNAM can be an integral part of a total reactor control system, where it can be used to continuously monitor and adjust the inlet core pressure so as to implement some predetermined sequence of reactivity insertion. The computer program DYNAM can also be used to modify system design (for example, reactor dimensions) and to observe the effect of such modifications on important system variables either during transient conditions, or during normal steady state operation. BPGCR Stability Analysis The BPGCR is a circulating fuel reactor using gaseous UFe as the fuel. A circulating fuel reactor by itself (i.e,, with no feedback effects) is an unstable configuration. This means that any perturbation acting on the reactor will make the power grow without limit as time progresses. It is the feedback effectsÂ— which are always present in a reactor operating at power, that make the system stable. Because of
PAGE 182
174 the gaseous nature of the fuel, plus the circulating fuel configuration, the BPGCR possesses some unique feedback effects which are found to render the system much more stable than conventional, solid core systems. To accurately predict the stability of the BPGCR, a nonlinear stability analysis of the system is performed using DYNAM, Some of the more important results obtained from this analysis are given below; a) the BPGCR is found to be an extremely stable configuration. Even after significant reactivity perturbations, the system returns to steady state within a short period of time. It is observed that after an external reactivity perturbation of the order of $1. the resulting power oscillations have amplitudes of no more than 30% of the steady state value. b) the inherent BPGCR response to external reactivity perturbations is found to be very different from that of a conventional, solid core, noncirculating fuel reactor. It is found that at full power, steady state conditions, the BPGCR power cannot be increased by merely applying an external positive reactivity step insertion. To increase the BPGCR power, it is necessary to provide additional positive reactivity during the transient by varying a system operating parameter (e.g., by increasing the core inlet pressure). c) variations in the external loop circulation time are found to influence the dynamic behavior of the BPGCR for the examined
PAGE 183
175 transients and the combination of core residence time and external loop circulation times considered in this dissertation. d) the most important feedback effect that influences the stability of the BPGCR is found to be the fuel mass (density) effect. Depending on the value of this fuel mass feedback coefficient, the BPGCR can be made stable, conditionally stable or unstable. e) limits on desired values of the fuel mass feedback coefficient are determined. It is found that for a certain range of values of this coefficient, the reactor exhibits oscillatory response to reactivity perturbations. But these oscillations are heavily damped and, by properly selecting the operating regime of the BPGCR, these oscillations can be fully eliminatedÂ— if such oscillations are deemed to be undesirable. If the BPGCR system is needed to provide pulsed mode power, controlled oscillations in the power output can help to reduce, or even eliminate, some of the power conditioning needs. In such a situation, the value of the fuel mass feedback coefficient can be carefully selected to provide sustained oscillatory output. Sometimes performing a nonlinear analysis can be tedious, and it is preferable to do a linear stability analysis if such an analysis can properly predict the stability of the system. The nonlinear system equations of the BPGCR are linearized and the stability of the BPGCR is determined using the established methods of linear stability analysis. Comparing the results of the linear analysis with that of the nonlinear analysis, it is established that for the most part, a linear
PAGE 184
176 stability analysis is able to qualitatively predict correctly the stability of the BPGCR. So, for initial system characterization, a linear stability analysis may prove useful. For the final design calculations, however, one needs to perform a fullfledged, nonlinear analysis to obtain accurate bounds on stability. Recommendations for Future Work Foundations for the dynamic analysis of the BPGCR are laid in this dissertation. Some refinementsÂ— which will make the analysis more accurate, are discussed below. Some of the recommendations are short term and others need a longer time to implement. On a short term basis, one can perform the foUoAving tasks. a) Some further parametric studies can be undertaken to study the effect of changing the system dimensions or design. For example, the inlet and exit throat areas can be varied, and the effect on system dynamics can be studied. These nozzle areas influence the fuel mass inflow and outflow rates from the core and, hence, will influence the reactivity variations during transients. The nozzles can also be replaced by valves or orifices to see whether they will make the BPGCR more stable or not. Similarly, some of the balanceofplant variables like radiator areas can be varied to see the influence of increasing or decreasing the external loop damping on the d)mamic behavior of the BPGCR. b) In order to reduce computational time, only one delayed neutron group is considered in DYNAM while solving for the
PAGE 185
177 neutronics equations. A sample case, comparing the BPGCR responses for the case Avith six delayed neutron groups versus that for the one delayed neutron group are given in Figs. 4.8 and 4.9. These figures indicate that, for the considered transient, the results from the two cases are almost identical and, hence, a one delayed neutron group treatment is adequate. For certain other transients, treatment of all six delayed neutron groups may be needed to accurately predict the reactor power behavior. So, for final design calculations, all transients should be run with the six delayed neutron groups option of DYNAM. c) In DYNAM, the external loop circulation time is kept fixed during a transient. It is not allowed to change, for example, as a function of the core exit velocity of the fuel. This can be modified by fixing the external loop length and letting the fuel exit velocity determine the external loop circulation time. This might influence the d 3 mamic response of the BPGCR by influencing the core inlet temperatures and delayed neutron precursors reintroduced into the core. d) Some other methods of controlling the BPGCR can also be investigated. In this dissertation, the BPGCR is controlled by means of the core inlet pressure. How the system will respond to control by means of the fuel mass flow rate into the core can also be investigated. Control of the BPGCR by means of core exit conditions (e.g., pressure) can also be studied, although this might be more difficult to
PAGE 186
178 implement in practice because of the extreme conditions (e.g., high fuel temperature) that exist at the core exit. e) In the balance ofplant, DYNAM models the MHD generator as a turbine generator. This can be refined by actually treating the MHD generator and seeing how the power oscillations inside the BPGCR affect the performance of the MHD generator and vice versa. f) It is found in this dissertation that variations in the fuel external loop circulation time influence the dynamic behavior (final steady state power, amplitude of initial power oscillations, etc.) of the BPGCR after an external reactivity perturbation. For a circulating fuel reactor, it is actually the ratio of the core residence time and the external loop circulation time that is important. This ratio decides what fraction of the delayed neutrons is lost from the system due to the circulation of delayed neutron precursors. So, other values for the ratio of these two times should also be examined. The way DYNAM is written currently, there is only provision for changing the external loop circulation time of the fuel (and not for changing the core residence time of the fuel). To study the effect of varying the core residence time of fuel on the BPGCR dynamics, DYNAM needs to be modified. One way to vary the core residence time of the fuel would be to artificially vary the length of the core and, hence, the transit time of the fuel through the core. Some of the refinements that can be incorporated into the BPGCR dynamic analysis on a long term basis are discussed next.
PAGE 187
179 a) In DYNAM, point reactor kinetics (PRK) equations are used to predict the behavior of system power during transients. PRK equations are capable of predicting the behavior of the system power with reasonable accuracy, but a better approach would be to use spacetime kinetics to predict the system power behavior. In DYNAM, a uniform gas density is assumed throughout the core. This is certainly not true. Unlike in a solid core system, the fuel density varies significantly within the core of the BPGCR. Since the fuel is flowing continuously through the core, the flow pattern and the flow rate of the fuel through the core significantly influences the temperature profile and, hence, the density profile of the fuel within the core. Another reason for the nonuniform density profile within the core is the fact that the thermal flux peaks near the beryllium moderator. This can result in an increased power production and, hence, increased heat generation near the walls of the BPGCR depending upon the fuel gas density distribution within the core. Note that an increased heat generation near the walls implies a lower gas density near the walls which, in turn, implies a reduced fission heat generation compared to the center of the core. In practice, this effect will be partly cancelled by any wall cooling that might be present, i.e., some coolant will be flowing through the walls to keep its temperature from exceeding the melting point of beryllium (or BeO). This will have the effect of making the fuel denser near the wall. It is apparent from the above discussion that the gas density
PAGE 188
180 distribution within the core can affect the spatial flux distribution within the core and vice versa. Another spacetime effect that might be present in the system is the variation of the core average flux relative to that of the reflector average flux during a transient. For a cavity configuration, like the BPGCR, the fission neutrons produced within the core get moderated in the reflector surrounding the core (no moderation takes place within the core proper). The mean time taken by the fast neutrons to return to the core, after slowing down in the reflector, is only of the order of a few milliseconds. But, when considering fast transients, the variation of the core average flux with respect to that of the reflector average flux needs to be considered. Such spacetime effects cannot be properly handled by PRK models. b) Another phenomenon that needs to be investigated is the effect of acoustic oscillations Avithin the BPGCR core. In a cavity configuration like that of the BPGCR, employing compressible fluids and with pressure dependent internal heat generation, it is an observed fact that acoustic oscillations are frequently generated. Such acoustic oscillations in the form of standing pressure waves have been predicted for gaseous core reactors. These acoustic oscillations can interact with the neutronic (amplitude) oscillations analyzed in this dissertaUon, or with any spacetime flux oscillations that might be present in the system (spacetime flux oscillations are normally
PAGE 189
181 observed only in looselycoupled cores). This can be a significant effect, and needs to be investigated. c) A third refinement in the analysis would be to include the second half of the BPGCR (in this dissertation, only one half of the BPGCR is considered), and introduce the coupling effects between the two halves while performing the dynamic analysis. Coupling between two such cores can be significant, and can affect the d 5 niamic response of both halves. Panicker (Ref.3) has shown that, for the pulsed gas core reactors (PGCRs), the coupling between these small cores can significantly influence their dynamic response.
PAGE 190
APPENDIX A THE COMPUTER PROGRAM DYNAM Introduction For analyzing the transition of the BPGCR from a low power mode of operation to a high power mode of operation, a computer program called DYNAM is developed in this dissertation. In Chapter 3, a detailed description of this program and the models used in it are given. DYNAM is also used to perform a nonlinear stability analysis of the BPGCR in Chapter 4. This Appendix contains a listing of the FORTRAN program. For a detailed description of the models used in this program, refer to Chapter 3 and Appendix B of this dissertation. 182
PAGE 191
183 PROGRAM DYNAM INTEGER N, KOUNT, K1 , KMAX, KOUNTl , KOUNT2 , NERR, KMAXl REAL*8 ATOMS (5) ,NUMDEN(4) ,DENSIT{5) , MASS (14) , H (7) , SECOND (8) , 1 KEFF , ENRICH, DEN ( 7 ) , POWER, DENEND ( 7 ) , MEANDN, CINLET ( 6 ) , 1 COUT (6,5500), MEANCP, MASSl , MASS2 , AVGMAS, MASSIN, MLEFT, 1 TOTMAS,NHE,NUF6, THROAT, VOL, HEIGHT, DIA, RADIUS, AREAIN, 1 AREAEX, AREAl , BETHIK, BEMASS, CPBE, RHOBE, RUNIV, RFUEL, 1 SBCONS,CPBYCV,EMISS,LSTAR,BETA(6) ,ALAMDA(6) ,MACHIN, 1 CORVEL, MACHSQ, INRATE, EXRATE, RHOTOT, RHCmX, GAMMA, 1 STAGPI, STAGPC, STAGTC, STAGTI, TINLET, TNEW, CORTMP, BETMP, 1 COREPR, RATIO, FACTOR, FACTSQ, CYLINl , CYLIN2 , TERM ( 3 ) , 1 DELTAT, HEATNG (3) , STEP, TLOOP, EXPO, TIME, TIME2,TIMAX,T2, 1 T3,S0UNDI,DENIN,CRITPR,H2(7) ,SECOND2(8) ,RHO(14) , 1 RHOEXT, RHOINT, TEX, TRAD, T01,T04,TEXIT (5500) ,RSLOPE, 1 EXITPR(5500) COMMON /INTPRM/ LSTAR, BETA, ALAMDA CCayiMON /CORDIA/ VOL, HEIGHT, DIA, RADIUS, AREAIN, AREAEX, THROAT, 1 BETHIK, ENRICH COMMON /CONST/ RUNIV, RFUEL, CPBYCV, MEANCP, CPBE, SBCONS, EMISS COMMON /TIMECO/ TIME, TIMAX, TLOOP, STEP C OPEN (UNIT=1,FILE='FINAL2.0P1',STATUS='NEW) C OPEN (UNIT=2,FILE=Â’FINAL2.0P2' ,STATUS=Â’NEW' ) C OPEN (UNIT=3,FILE='FINAL2.0P3' ,STATUS='NEW' ) OPEN (UNIT=4,FILE='FINAL2.0P4',STATUS='NEW ) DATA ATOMS /I . 2712, 0 . 2243, 8 . 9730, 19 . 8688, 0 . 1214/ DATA MASS /1 . 0, 2 . 0, 3 . 0, 4 , 0, 5 . 0, 6 . 0, 7 . 0, 8 , 0, 9 . 0, 9 . 2, 9 . 4, 1 9.6,9.8,10.0/ DATA KOUNT, KOUNT2,Kl, KMAX, N /O, 0, 0, 1000, 7/ DATA DEN (1 ), CINLET /O . 90E+8, 6*0 . 0/ DATA STAGPI, STAGTI, CORTMP, TEX, TRAD, TOl, BETMP, CORVEL /35.0, 1 1200.0,1500.0,2000.0,1200.0, 900.0, 800.0,30.0/ DATA RHOTOT, GAMMA, RHOMAX, TIME2/0 . 003, 0 . 02, 0 . 00710,0 . 0/ RHOBE = 1800.0 RSLOPE =0.05
PAGE 192
10 20 25 30 50 184 DO 10 1=1, 4,1 NUMDEN(I) = ATCa^Sd) / (VOL * l.OE+6) CONTINUE DENSIT(5) = AT01S(5) DO 25 J=l,14,l DO 20 I = 1,4,1 DENSIT(I) = NUMDEN(I) * MASS ( J) CONTINUE CALL CRITIC (DENSIT,KEFF) RHO(J) = 1.0 (1.0 / KEFF) CONTINUE CALL TABLE (RHO,MASS,RHOTOT,MASSl) DO 30 I = 2,N, 1 DEN (I) = BETA(Il) *DEN(1) / (LSTAR*ALAMDA(I1) * (l.ORHOTOT) ) CONTINUE CYLINl = DIA / 2.0 CYLIN2 = CYLINl + (BETHIK / 100.0) BEMASS = RHOBE * 3.1416 * (CYLIN2 ** 2 CYLINl **2) * HEIGHT AREAl = 3.1416 * DIA * HEIGHT CONTINUE TIME = TIME + STEP IF (TIME .GE. 1.5) THEN GO TO 65 ELSE IF (TIME .GE. TLOOP) THEN RHOTOT = 0.0035 ENDIF CALL RK4 (DEN, CINLET,N, RHOTOT, CORVEL,DENEND) MEANDN = (DEN(l) + DENEND(l)) / 2.0 IF (DENEND(l) .GT. l.OE+9) THEN PRINT *, 'POWER EXCEEDED 1000 MWS' GO TO 90
PAGE 193
185 ENDIF NUF6 = MASSl * 0.868 / 0.3495 NHE = MASSl * 0.132 / 0.004 COREPR = (NUF6 + NHE) * RUNIV * CORTMP * 9.872E6 / VOL CALL ACTUAL (STAGPI, STAGTI, COREPR, CPBYCV, TINLET,MACHIN, NERR) IF (NERR .GT. 0) THEN KOUNT2 = KOUNT2 + 1 PRINT*, 'AT', TIME, 'SECS', KOUNT2, ' CORE INLET 1 STAGNATION PRESSURE (S) LOWER THAN COREPR! ! ! ! ' ENDIF CORVEL = 0.1 * MACHIN * DSQRT (CPBYCV * RFUEL * TINLET) STAGPC = STAGPI C CALL FLORTE (CORTMP, COREPR, EXRATE) EXRATE = (0.6847 * 1.013E+5 * COREPR * THROAT) / (DSQRT (RFUEL 1 * CORTMP) ) C CALL REACT (TIME2, ALAMDA, STEP, GAMMA, RHOMAX, RHOTOT) CALL TABLE (RHO, MASS, RHOTOT, MASS2) INRATE = EXRATE + ( (MASS2 MASSl) / STEP) CALL CONTRL (INRATE, COREPR, TINLET, AREAIN, STAGPI) CALL TEMP ( TINLET , INRATE , EXRATE , STEP , AREAl , BEMASS , DENEND ( 1 ) , 1 MASS2, CORVEL, COREPR, TEX, CORTMP, BETMP, HEATNG) IF (TIME .GE. TLOOP) THEN DO 52 I = 1,6,1 CINLET(I) = COUT(I,l) * DEXP (ALAMDA (I ) * TLOOP) 52 CONTINUE CALL INTEMP (TEXIT (1) , T04, TOl, TRAD, STAGTI, DELTAT)
PAGE 194
186 IF (STAGTI .LE. 600.0) THEN STAGTI = 600.0 ENDIF DO 54 I = 1,6,1 DO 53 J = 1, (KL1) ,1 GOUT (I, J) = COUT(I,J+l) 53 CONTINUE 54 CONTINUE DO 55 I = 1, (Kl1) ,1 TEXIT(I) = TEXIT(I+1) 55 CONTINUE DO 56 I = 1,6,1 COUT(I,Kl) = DENEND(I+1) 56 CONTINUE TEXIT(Kl) = TEX ELSE K1 = K1 + 1 DO 57 I = 1,6,1 COUT(I,Kl) = DENEND(I+1) 57 CONTINUE TEXIT(Kl) = TEX ENDIF KOUNT = KOUNT + 1 IF (KOUNT .GE. KMAX) THEN C WRITE (1,92) TIME,RHOTOT,MASS2,DENEND (1) C WRITE (2, 94) TIME, COREPR, STAGPI, CORTMP, STAGTI, CORVEL, EXRATE C WRITE (3, 96) TIME, TEX, T04, TOl, TRAD, STAGTI, DELTAT KOUNT = 0 ENDIF DO 60 I = 1,N DEN (I) = DENEND(I)
PAGE 195
187 60 CONTINUE MASSl = MASS2 GO TO 50 65 CONTINUE TIME =0.0 C DO 68 I = 1,K1,1 C EXITPR(I) = STAGPI C 68 CONTINUE C STAGPI = STAGPI +1.0 RHOINT = RHOTOT C WRITE(1,97) TIME, RHOTOT, RHOINT,MASS2,DENEND(l) WRITE(4, 98) TIME, TIME, RHOINT, DENEND (1) , STAGPI C WRITE (2, 94) TIME, COREPR, STAGPI, CORTMP, STAGTI, CORVEL, EXRATE KMAX = 200 KOUNT = 0 KMAXl = 500 KOUNTl = 0 RHOEXT = 0.00058 RHOTOT = RHOTOT + RHOEXT 70 CONTINUE TIME = TIME + STEP IF (TIME .GT. TIMAX) THEN GO TO 90 ENDIF CALL RK4 (DEN, CINLET, N, RHOTOT, CORVEL, DENEND ) C MEANDN = (DEN(l) + DENEND(l)) / 2.0 C IF (DENEND (1) .GT. l.OE+9) THEN C PRINT *, 'POWER EXCEEDED 1000 MWS'
PAGE 196
188 C GO TO 90 C ENDIF NUF6 = MASSl * 0.868 / 0.3495 NHE = MASSl * 0.132 / 0.004 COREPR = (NUF6 + NHE) * RUNIV * CORTMP * 9.872E6 / VOL CALL ACTUAL (STAGPI , STAGTI , COREPR, CPBYCV, TINLET, MACHIN, NERR) IF (NERR .GT. 0) THEN KOUNT2 = KOUNT2 + 1 PRINT*, 'AT', TIME, 'SECS', KOUNT2, ' CORE INLET 1 STAGNATION PRESSURE (S) LOWER THAN COREPR! ! ! ! ' ENDIF CORVEL = 0.1 * MACHIN * DSQRT (CPBYCV * RFUEL * TINLET) STAGPC = STAGPI C CALL FLORTE (MACHSQ, FACTSQ, CORTMP, STAGPC, STAGTC, EXRATE) EXRATE = (0.6847 * 1.013E+5 * COREPR * THROAT) / (DSQRT 1 (RFUEL * CORTMP) ) SOUNDI = DSQRT (CPBYCV * RFUEL * TINLET) DENIN = COREPR * 1.013E+5 / (RFUEL * TINLET) VELIN = SOUNDI * MACHIN CRITPR = 0.528 * STAGPI IF (COREPR .LE. CRITPR) THEN INRATE = (0.6847 * 1.013E+5 * STAGPI * AREAIN) / 1 DSQRT (RFUEL * STAGTI) ELSE INRATE = DENIN * AREAIN * VELIN ENDIF MASS2 = MASSl + STEP * (INRATE EXRATE) CALL TABLE (MASS, RHO, MASS2, RHOINT) C RHOINT = RHOINT + RSLOPE * (MASS2 MASSl) CALL TEMP (TINLET, INRATE, EXRATE, STEP, AREAl, BEMASS, DENEND (1 ) , MASS2, CORVEL, COREPR, TEX, CORTMP, BETMP, HEATNG) 1
PAGE 197
189 DO 72 I = 1,6,1 CINLET(I) = GOUT (1,1) * DEXP ( ALAMDA ( I ) * TLOOP) 72 CONTINUE CALL INTEMP (TEXIT (1 ) , T04, TOl, TRAD, STAGTI, DELTAT) DO 74 I = 1,6,1 DO 73 J = 1, (Kll),l COUT(I,J) = COUT(I,J+l) 73 CONTINUE 74 CONTINUE DO 75 I = 1, (Kl1) ,1 TEXIT (I) = TEXIT (I+l) C EXITPR(I)= EXITPR(I+1) 75 CONTINUE DO 76 I = 1,6,1 COUT(I,Kl) = DENEND(I+1) 76 CONTINUE TEXIT (Kl) = TEX C EXITPR(Kl) = COREPR C KOUNT = KOUNT + 1 C IF (KOUNT .GE. KMAX) THEN C WRITE(1,97) TIME,RHOTOT,RHOINT,MASS2,DENEND (1) C WRITE (2, 94) TIME, COREPR, STAGP I, CORTMP, STAGTI, CORVEL, EXRATE C WRITE (3, 96) TIME, TEX, T04, TOl, TRAD, STAGTI, DELTAT C KOUNT = 0 C ENDIF KOUNTl = KOUNTl + 1 IF (KOUNTl .GE. KMAXl) THEN WRITE (4, 98) TIME, TIME, RHOINT, DENEND (1) , STAGPI KOUNTl = 0 ENDIF DO 80 I = 1,N
PAGE 198
190 DEN (I) = DENEND(I) 80 CONTINUE MASSl = MASS2 RHOTOT = RHOINT + RHOEXT C STAGPI = EXITPR(l) GO TO 70 90 CONTINUE 92 FORMAT(2F10.5,F10.3,2E12.4) 94 FORMAT(F10.4,6F10.3) 96 FORMAT(F10.4, 6F10.3) 97 FORMAT(3F10.5,F10.3,E12.5) 98 FORMAT(F10.3,1X,1H, ,F10.5,1X,1H, ,F10.5,1X,1H, ,E12.3,1X,1H, ,F8.2) END SUBROUTINE CRITIC (DENSIT, KEFF) C THIS SUBROUTINE CALCULATES THE CORE EFFECTIVE MULTIPLICATION C FACTOR FOR A GIVEN NUMBER DENSITY OF FUEL AND BERYLLIUM REFLECTOR. C IT USES A 2 GROUP, 2 REGION DIFFUSION THEORY CALCULATION ON A C SPHERICAL SYSTEM FOR ESTIMATING THE SYSTEM KEFF. THE MICROSCOPIC C GROUP CONSTANTS WHICH THE PROGRAM USES FOR CALCULATING THE SYSTEM C KEFF ARE OBTAINED VIA INDEPENDENT XSDRNPM CALCULATIONS AND ARE FED C TO THE PROGRAM AS INPUT. SO, IF THERE ARE ANY DRASTIC CHANGES IN THE C SYSTEM GECMETRY OR FUEL CONCENTRATION, THESE CONSTANTS NEED TO BE C REEVALUATED . INTEGER INTRVL(2) ,N1,N2,N REAL*8 SORCEl (400) ,SORCE2 (400) ,GP2SRC(400) ,GP1FLX(400) , 1 GP2FLX (400) , DENSIT (5) ,TOTALl (5) ,TOTAL2 (5) ,ABS1(5) , 1 ABS2 (5) ,FISS1 (5) ,FISS2 (5) , TRANSl (5) ,TRANS2 (5) ,ONET01 (5) , 1 ONET02 (5) ,TWOT02 (5) , SIGFl (2) , SIGF2 (2) , SIGRl (2) , SIGR2 (2) , 1 SIGAl (2) , SIGA2 (2) , SIGH (2) , SIG12 (2) , SIG22 (2) , SIGTRl (2) , 1 SIGTR2,D1 (2) ,D2 (2) , STEP (2) , NUl (2) ,NU2 (2) , NU (2) , ENRICH, ^ KOLD, KNEW, RADIUS, BETHIK,BTHIK2, EPS, TESTl,OLDSUM, 1 NEWSUM,TERM(6) ,A1 (400) ,B1 (400) , Cl (400) , El (400) , FI (400) ,
PAGE 199
1 1 A2 (400) ,B2 (400) ,C2 (400) , E2 (400) ,F2 (400) , KEFF, VOL, HEIGHT DIA, AREAIN, AREAEX, THROAT 191 1 1 1 CCMION /MICRO/ T0TAL1,T0TAL2,FISS1,FISS2,ABS1,ABS2,NU1,NU2, TRANSl , TRANS2 , ONETOl , ONET02 , TWOT02 COMMON /MACRO/ Dl, D2, SIGAl, SIGA2, SIGFl, SIGF2, SIGH, SIG12, SIG22, SIGTRl, SIGTR2, SIGRl, SIGR2, NU COMMON /CORDIA/ VOL, HEIGHT, DIA, RADIUS, AREAIN, AREAEX, THROAT, BETHIK, ENRICH DATA INTRVL,KOLD,EPS /70, 200, 1 . 0, 1 . OE4/ CALL GRPCON (DENSIT, ENRICH) BTHIK2 = BETHIK + (0.71 / SIGTR2) STEP(l) = RADIUS / INTRVL(l) STEP (2) = BTHIK2 / INTRVL(2) N1 = INTRVL(l) N2 = INTRVL(2) N = N1 + N2 CALL COEFF(N,Nl, STEP, RADIUS, Dl (1) ,D2 (1) , SIGRl (1) ,SIGR2 (1) ,A1,B1, 1 Cl, El, FI) CALL COEFF (N,N1, STEP, RADIUS, Dl (2) ,D2 (2) , SIGRl (2) ,SIGR2(2) ,A2,B2, 1 C2,E2,F2) DO 110 I = 1,N1,1 SORCEl(I) =1.0 110 CONTINUE DO 120 I = (Nl+1) , (N1) ,1 SORCEl(I) =0.0 120 CONTINUE TERM(l) = (4.0 /3.0) * 3.1416 * ( STEP(l) / 2.0)**3 TERM(2) = (3.1416/3.0) * (STEP(1)**3) TERM(3) = 3.1416 * ((2.0 * STEP(l) * RADIUS * RADIUS) (RADIUS * STEP(l) *STEP(D) + (STEP(1)**3 / 6.0)) 1
PAGE 200
OLDSUM = SORCEl(l) * TERM(l) DO 130 1=1, (Nll),l OLDSUM = OLDSUM + SORCEl(I) * TERM(2) * (1,0 + 12.0 * I * 130 CONTINUE OLDSUM = OLDSUM + SORCEl(Nl) * TERM (3) 150 CONTINUE DO 155 I = 1, (N1) ,1 SORCEl(I) = SORCEl(I) / KOLD 155 CONTINUE CALL GRPFLX (N, N1 , SORCEl , A1 , B1 , Cl , El , FI , GPIFLX) DO 160 1=1, (Nll),l GP2SRC(I) = SIG12(1) * GPIFLX (I) 160 CONTINUE GP2SRC(N1) = (SIG12(1) + SIG12(2)) / 2.0 * GPIFLX (Nl) DO 170 I = (Nl+1), (Nl),l GP2SRC(I) = SIG12(2) * GPIFLX (I) 170 CONTINUE CALL GRPFLX (N, Nl , GP2SRC, A2 , B2 , C2 , E2 , F2 , GP2FLX ) TERM(4) = NU(1) * SIGFl(l) TERM(5) = NU(2) * SIGFl (2) DO 175 I = 1,N1,1 SORCE2(I) = TERM(4) * GPIFLX(I) + TERM(5) * GP2FLX(I) 175 CONTINUE DO 180 I = (Nl+1) , (N1) ,1 SORCE2(I) =0.0 180 CONTINUE NEWSUM = SORCE2(l) * TERM(l) DO 185 I = 1, (Nl1) ,1 NEWSUM = NEWSUM + SORCE2(I) * TERM(2)*(1.0 + 12.0 * I * I) 185 CONTINUE
PAGE 201
193 NEWSUM = NEWSUM + (S0RCE2(N1) * TERM(3) ) KNEW = KOLD * NEWSUM / OLD SUM TESTl = DABS ({KNEW KOLD) / KNEW) IF (TESTl .LE. EPS ) THEN GO TO 200 ENDIF DO 190 I = 1, (N1) ,1 SORCEl(I) = SORCE2(I) 190 CONTINUE KOLD = KNEW OLDSUM = NEWSUM GO TO 150 200 CONTINUE KEFF = KNEW RETURN END SUBROUTINE GRPCON (DENSIT, ENRICH) C THIS SUBROUTINE CALCULATES THE MACROSCOPIC GROUP CONSTANTS WHICH C ARE NEEDED BY THE CRITICALITY PROGRAM FOR ESTIMATING THE SYSTEM C KEFF. REAL*8 ENRICH, DENSIT (5) ,TOTALl (5) ,TOTAL2 (5) ,ABS1 (5) ,ABS2(5) , 1 TRANSl (5) ,TRANS2 (5) , FISSl (5) , FISS2 (5) ,ONET01 (5) , 1 ONET02 (5) ,TWOT02 (5) , SIGTOT (2) , SIGRl (2) ,D1 (2) ,D2 (2) , 1 SIGR2 (2) , SIGAl (2) , SIGA2 (2) , SIGFl (2) , SIGF2 (2) , SIGTRl (2) , 1 SIGTR2,SIG11 (2) ,SIG12(2) ,SIG22 (2) ,NU1 (2) ,NU2 (2) ,NU(2) COMMON /MICRO/ TOTALl, TOTAL2, FISSl, FISS2, ABS1,ABS2, NUl, NU2, 1 TRANSl , TRANS2 , ONETOl , ONET02 , TWOT02 COMMON /MACRO/ D1,D2, SIGAl, SIGA2, SIGFl, SIGF2, SIGH, SIG12, SIG22, 1 SIGTRl, SIGTR2, SIGRl, SIGR2,NU DO 200 I = 1,2 SIGAl (I) =0.0
PAGE 202
194 SIGFl(I) =0.0 SIGTR1(I)= 0.0 200 CONTINUE SIGll(l) = 0.0 SIG12(1) = 0.0 SIG22(1) =0.0 SIGTOT(l) =0.0 DO 250 I = 1,4,1 SIGTRl(l) = SIGTRl(l) + DENSIT(I) * TRANSl(I) SIGTR1(2) = SIGTRK2) + DENSIT(I) * TRANS2(I) SIGAl(l) = SIGAl(l) + DENSIT(I) * ABSl(I) SIGA1(2) = SIGA1(2) + DENSIT(I) * ABS2(I) SIGFl(l) = SIGFl(l) + DENSIT(I) * FISSl(I) SIGF1(2) = SIGF1(2) + DENSIT(I) * FISS2(I) SIGll(l) = SIGll(l) + DENSIT(I) * ONETOl (I) SIG12(1) = SIG12(1) + DENSIT(I) * ONET02(I) SIG22(1) = SIG22(1) + DENSIT(I) * TWOT02(I) SIGTOT(l) = SIGTOT(l) + DENSIT(I) * TOTALl(I) 250 CONTINUE Dl(l) = 1.0 / (3.0 * SIGTRl(D) Dl(2) = 1.0 / (3.0 * SIGTR1(2)) D2(l) = 1.0 / (3.0 * DENSIT(5) * TRANSl (5) ) D2(2) = 1.0 / (3.0 * DENSIT(5) * TRANS2(5)) SIGA2(1) = DENSIT(5) * ABSl (5) SIGA2(2) = DENSIT(5) * ABS2(5) SIGF2(1) =0.0 SIGF2(2) =0.0 SIGTR2 = DENSIT(5) * (TRANSl (5) + TRANS2 (5) ) / 2.0 SIGH (2) = DENSIT(5) * ONETOl (5) SIG12(2) = DENSIT(5) * ONET02 (5) SIG22(2) = DENSIT(5) * TW0T02(5) SIGTOT(2) = DENSIT(5) * TOTALl (5) SIGRl(l) = SIGTOT(l) SIGll(l) SIGRK2) = SIGA1(2) SIGR2(1) = SIGTOT(2) SIGH (2) SIGR2(2) = SIGA2(2)
PAGE 203
195 NU(1) = ENRICH * NUl(l) + (1.0 ENRICH) * NU1(2) NU(2) = ENRICH * NU2(1) + (1.0 ENRICH) * NU2(2) RETURN END SUBROUTINE COEFF (N,N1, STEP, RADIUS, Dl, D2, SIGl, SIG2, A, B, C, E, F) C THIS SUBROUTINE CALCULATES THE ELEMENTS OF THE COEFFICIENT MATRIX C WHICH IS USED FOR SOLVING THE GROUP FLUXES (BOTH FAST AND THERMAL) . C THIS INFACT RESULTS IN A PENTA DIAGONAL COEFFICIENT MATRIX. INTEGER N,N1 REAL*8 STEP(2),A(400) ,B(400) ,C(400) ,E(400) ,F(400) ,D1,D2,SIG1, 1 SIG2, DENOM, FACT (7) , RADIUS, RI DO 300 I = 1,N A(I) = 0.0 B(I) = 0.0 C(I) = 0.0 E(I) = 0.0 F(I) = 0.0 300 CONTINUE FACT(l) = Dl / (2.0 * STEP(l)) FACT (2) = D2 / (2.0 * STEP (2)) A(N1) = FACT(l) B(N1) = 4.0 * FACT(l) C(N1) =3.0 * (FACT(l) + FACT (2)) E(N1) = 4.0 * FACT (2) F(N1) = FACT (2) FACT (3) = Dl / STEP(1)**2 C(l) =FACT(3) * 2.0 SIGl E(l) = FACT (3) * 2.0 DO 320 1=2, (Nl1) ,1 DENOM =1.0/1
PAGE 204
B(I) = FACTO) * (1.0 DENOM) E(I) = FACTO) * (1.0 + DENOM) C(I) =FACT(3) * 2.0 SIGl 320 CONTINUE FACT (4) = D2 / STEP (2) DO 340 I = (Nl+1) , (N1) ,1 RI = RADIUS + (I Nl) * STEP (2) FACT (5) = (1.0 / STEP (2) 1.0 / RI) FACT (6) = (1.0 / STEP (2) + 1.0 / RI) B(I) = FACT (4) * FACT (5) E(I) = FACT (4) * FACT (6) C(I) =2.0 * FACT (4) / STEP (2) SIG2 340 CONTINUE DO 360 1=2, (N2) ,1 B(I) = B(I) / C(Il) A(I+1) = A(I+1) / C(Il) C(I) = C(I) B(I) * E(Il) E(I) = E(I) B(I) * F(Il) B(I+1) = B(I+1) A(I+1) * E(Il) C(I+1) = C(I+1) A(I+1) * F(Il) 360 CONTINUE B(Nl) = B(Nl) / C(N2) C(Nl) = C(Nl) B(Nl) * E(N2) RETURN END SUBROUTINE GRPFLX (N, Nl, R, A, B, C, E, F, X) C THIS SUBROUTINE SOLVES FOR THE GROUP FLUXES USING THE COEFFICIENT C MATRIX CALCULATED BY THE PREVIOUS SUBROUTINE. THE SOLUTION C TEC NIQUE IS VERY SIMILAR TO THE THOMAS ALGORITHM WHICH IS USED C FOR SOLVING TRIDIAGONAL MATRICES. IN OUR CASE BECAUSE OF THE C SPECIAL INTERFACE CONDITON, WE END UP WITH A PENTAD lAGONAL MATRIX.
PAGE 205
197 INTEGER N,N1 REAL*8 A (N1 ) , B (N1 ) , C (N1 ) , E (N1 ) , F (N1 ) , R (N1 ) , X (N1 ) DO 400 1 = 1, (N1) ,1 R(I) = R(I) 400 CONTINUE R(N1) = 0.0 DO 420 I = 2, (N2) ,1 R(I) = R(I) B(I) * R(Il) R(I+1) = R(I+1) A(I+1) * R(Il) 420 CONTINUE R(Nl) = R(Nl) B(Nl) * R(N2) X(Nl) = R(Nl) / C(Nl) X(N2) =(R(N2) E(N2) * X(Nl)) / C(N2) DO 450 I = (N3),l,l X(I) = (R(I) X(I+1) * E(I) X(I+2) * F(D) / C(I) 450 CONTINUE RETURN END SUBROUTINE SPLINE (RHO, MASS, H, SECOND) C THIS SUBROUTINE SOLVES FOR THE SECOND DERIVATIVES NEEDED FOR CALCUC DATING THE COEFFICIENTS OF A CUBIC SPLINE FIT. REAL*8 RH0(8),MASS(8),H(7),TERM(7),A(6),B(6),C(6),R(6),G{6), 1 SECOND (8) DO 500 I = 1,7,1 H(I) = RHO(I+l) RHO(I) TERM(I) = (MASS(I+1) MASS (I)) / H(I) 500 CONTINUE DO 520 I = 1,6,1 A(I) = H(I)
PAGE 206
198 B(I) = 2.0 * ( H(I+1) + H(I) ) C(I) = H(I+1) R(I) =(TERM(I+1) TERM(D) * 6.0 520 CONTINUE A(l) = 0.0 C(6) = 0.0 DO 540 I = 2,6,1 A(I) = A(I) / B(Il) B(I) = B(I) A(I) * C(Il) R(I) = R(I) A(I) * R(Il) 540 CONTINUE G(6) = R(6) / B(6) DO 550 I = 5, 1,1 G(I) = (R(I) G(I+1) * C(I) ) / B(I) 550 CONTINUE SECOND (1) =0.0 SECOND (8) =0.0 DO 580 I = 1,6,1 SECOND (I+l) = G(I) 580 CONTINUE RETURN END SUBROUTINE TABLE (X, Y, IN, OUT) C THIS SUBROUTINE CALCULATES THE FUEL MASS IN THE CORE CORRESPONDING C TO A GIVEN REACTIVITY OR THE REACTIVITY CORRESPONDING TO A GIVEN C CORE FUEL MASS BY LOOKING UP A TABLE. A LINEAR INTERPOLATION IS USED C FOR INTERPOLATING BETWEEN THE TABULATED VALUES. INTEGER I REAL*8 X(14) ,Y(14) ,IN,OUT IF (IN .LE. X(2)) THEN 1 = 1 ELSE IF (IN .LE. X(3)) THEN
PAGE 207
199 I = 2 ELSE IF (IN .LE. X(4)) THEN I = 3 ELSE IF (IN .LE. X(5)) THEN I = 4 ELSE IF (IN .LE. X(6)) THEN I = 5 ELSE IF (IN .LE. X(7)) THEN I = 6 ELSE IF (IN .LE. X(8)) THEN I = 7 ELSE IF (IN .LE. X(9)) THEN I = 8 ELSE IF (IN .LE. X(10)) THEN I = 9 ELSE IF (IN .LE. X(ll)) THEN I = 10 ELSE IF (IN .LE. X(12)) THEN I = 11 ELSE IF (IN .LE. X(13)) THEN I = 12 ELSE I = 13 ENDIF OUT = Y(I) + (((Y(I+1) Y(I) )/(X(I+l) X(I))) * (INX(I))) RETURN END SUBROUTINE RK4 (BENIN, CIN, N, RHOTOT,VEL, DENOUT) C THIS SUBROUTINE PERFORMS FOURTH ORDER RUNGEKUTTA INTEGRATION WITH C CONSTANT TIME STEPS. INTEGER N REAL*8 DENIN(7) ,CIN(6) ,DENOUT(7) ,TEMP(7) ,SLOPEl (7) ,SLOPE2 (7) , 1 SLOPES (7) ,SLOPE4 (7) ,VEL, TIME, TIMAX, TLOOP, STEP, HFSTEP, 1 LSTAR,BETA(6) ,ALAMDA(6) ,RHOTOT, INVKEF,GENTME COMMON /INTPRM/ LSTAR, BETA, ALAMDA
PAGE 208
COMMON /TIMECO/ TIME, TIMAX, TLOOP, STEP INVKEF = 1.0 RHOTOT GENTME = LSTAR * INVKEF HFSTEP = STEP * 0.5 CALL RHS (BENIN, CIN,N, RHOTOT, GENTME, VEL, SLOPE!) DO 600 1=1, N TEMP(I)=DENIN(I) + HFSTEP*SLOPEl (I) 600 CONTINUE CALL RHS (TEMP, CIN, N, RHOTOT, GENTME, VEL, SLOPE2 ) DO 620 1=1, N TEMP(I)=DENIN(I) + HFSTEP*SLOPE2 (I ) 620 CONTINUE CALL RHS (TEMP, CIN, N, RHOTOT, GENTME, VEL, SLOPES) DO 640 1=1, N TEMP(I)=DENIN(I) + STEP*SLOPE3 (I) 640 CONTINUE CALL RHS ( TEMP , CIN, N, RHOTOT , GENTME , VEL , SLOPE 4 ) DO 650 1=1, N DENOUT(I)=DENIN(I)+(STEP/6.)*(SLOPEl(I) + 2.* 1 SLOPE2(I) + 2.*SLOPE3(I) + SLOPE4(I)) 650 CONTINUE RETURN END SUBROUTINE RK42 (IN,T04) C THIS SUBROUTINE PERFORMS FOURTH ORDER RUNGEKUTTA INTEGRATION WITH C CONSTANT TIME STEPS. INTEGER N REAL*8 IN(2) ,T04,TEMP(2) , SLOPE 1 (2) , SLOPE2 (2) , 1 SLOPES (2) , SLOPE 4 (2) , TIME, TIMAX, TLOOP, STEP, HFSTEP COMMON /TIMECO/ TIME, TIMAX, TLOOP, STEP N = 2 HFSTEP = STEP * 0.5 CALL RHS2 (IN, N, TO 4, SLOPE!)
PAGE 209
201 DO 660 I = 1,N TEMP (I) = IN(I) + HFSTEP * SLOPEl(I) 660 CONTINUE CALL RHS2 (TEMP,N,T04, SLOPE2) DO 670 I = 1,N TEMP (I) = IN(I) + HFSTEP * SLOPE2(I) 670 CONTINUE CALL RHS2 (TEMP, N, TO 4, SLOPES) DO 680 I = 1,N TEMP (I) = IN(I) + HFSTEP * SLOPES (I) 680 CONTINUE CALL RHS2 (TEMP, N, T04, SLOPE4 ) DO 690 I = 1,N IN(I) = IN(I) + (STEP/6.0) * (SLOPE! (I) + 2.0 * 1 SLOPE2(I) + 2.0 * SLOPES (I) + SLOPE4(I)) 690 CONTINUE RETURN END SUBROUTINE RHS (DENIN, CIN, N, RHO, GENTME, VEL, DNDT) C THIS SUBPROGRAM COMPUTES THE R.H.S OF THE SET OF NEUTRON AND C PRECURSOR DENSITY EQUATIONS. INTEGER N REAL*8 DENIN(N) ,CIN(6) ,DNDT(N) ,LSTAR,BETA(6) ,ALAMDA(6) , 1 RHO, GENTME, VEL, TIME, TIMAX, TLOOP, STEP, FLOTRM(6) , SUM COMMON /INTPRM/ LSTAR, BETA, ALAMDA COMMON /TIMECO/ TIME, TIMAX, TLOOP, STEP SUM =0.0 DO 700 I = 1,6,1 FLOTRM(I) = (VEL / 2.0) * (CIN(I) DENIN(I+1)) 700 CONTINUE DO 720 I = 1,6,1 SUM = SUM + (ALAMDA (I) * DENIN (I+l))
PAGE 210
202 720 CONTINUE DO 750 1=1, N,1 IF (I .EQ. 1) THEN DNDT(I) = ( (RHO0.0065) * DENIN(I)) / GENTME + SUM ELSE DNDT(I) = (BETA(Il) * DENIN(l) / GENTME) (DENIN(I) * 1 ALAMDA(ID) + FLOTRM(Il) ENDIF 750 CONTINUE RETURN END SUBROUTINE RHS2 (IN, N, T04, DNDT) C THIS SUBPROGRAM COMPUTES THE R.H.S OF THE SET OF HEAT EXCHANGER EXIT C TEMPERATURE AND AVERAGE RADIATOR TEMPERATURE EQUATIONS. INTEGER N, I REAL*8 IN(2) ,DNDT(2) ,T04 DO 800 1=1, N IF (I .EQ. 1) THEN DNDT(I) = 2.25 * IN(I) + 2.50 * IN(I+1) 0.250 * T04 ELSE DNDT(I) = (0.3676 * IN(Il)) (5.88E11 * (IN(I)**4)) 1 (0.7353 * IN(D) + (0.3676 * T04) ENDIF 800 CONTINUE RETURN END C SUBROUTINE REACT (TIME2, ALAMDA, STEP, SLOPE, RHOMAX, RHO) C C THIS SUBROUTINE CALCULATES THE TOTAL REACTIVITY IN THE CORE FOR
PAGE 211
203 C EACH TIME STEP C C REAL*8 TIME2, EXPO, ALAMDA,RH0,RH01,RH02,RH0MAX, STEP, SLOPE C C IF (TIME2 .GT. 0.0) THEN C EXPO = DEXP ((ALAMDA * TIME2) ) C RHOl = RHOMAX * EXPO C RH02 = 0.006458 * (1.0 EXPO) C RHO = RHOl + RH02 C TIME2 = TIME2 + STEP C ELSE C IF (RHO .GE. 0.0) THEN C SLOPE = 0.00015 C ENDIF C RHO = RHO + (SLOPE * STEP) C IF (RHO .GE. RHOMAX) THEN C TIME2 = TIME2 + STEP C ENDIF C ENDIF C C RETURN C END C c c c c c c c c c c c c SUBROUTINE FLORTE (MSQ, MACHEX, TEMP, STAGPC, STAGTC, MEX) THIS SUBROUTINE CALCULATES THE FUEL MASS OUTFLOW RATE FOR A GIVEN CORE STAGNATION CONDITIONS AND BACK PRESSURE. REAL*8 MSQ, MACHEX, MEX, PBACK, TEMP, STAGTC, STAGPC, 1 TEXIT, VOL, AREAIN, AREAEX, THROAT, RUNIV, RFUEL, CPBYCV, 1 BETHIK, ENRICH, HEIGHT, DIA, RADIUS, MEANCP, CPBE, RATIOl , 1 RATIO, SBCONS, EMISS COMMON /CONST/ RUNIV, RFUEL, CPBYCV, MEANCP, CPBE, SBCONS, EMISS COMMON /CORDIA/ VOL, HEIGHT, DIA, RADIUS, AREAIN, AREAEX, THROAT, BETHIK, ENRICH 1
PAGE 212
204 c c c MEX = (0.6847 * 1.013E+5 * STAGPC * THROAT) / (DSQRT (RFUEL * STAGTC) ) C 1 c RETURN C END SUBROUTINE ACTUAL (PO, TO, PI, CPBYCV, T1,MACH,NERR0R) C THIS SUBROUTINE CALCULATES THE ACTUAL CONDITIONS FOR A GIVEN C STAGNATION CONDITIONS . INTEGER NERROR REAL* 8 XI , X2 , RATIO, MSQ, MACH, PO , TO , PI , T1 , CPBYCV XI = CPBYCV 1.0 X2 = XI / CPBYCV RATIO = PO / PI MSQ = ((RATIO**X2) 1.0) * 2.0 / XI IF (MSQ .LT. 0.0) THEN MSQ = 1.0 * MSQ NERROR = 1 ELSE NERROR = 0 ENDIF MACH = DSQRT (MSQ) T1 = TO / (1.0 + (XI * MSQ / 2.0)) RETURN END SUBROUTINE STAG (M, T1 , CPBYCV, TO ) C THIS SUBROUTINE CALCULATES THE STAGNATION TEMPERATURE GIVEN THE C ACTUAL CONDITIONS AND THE SQUARE OF THE MACH NUMBER. REAL*8 M,T1, TO, CPBYCV, Y1,Y2, EXPO Y1 = CPBYCV 1.0 Y2 = Y1 * M / 2.0
PAGE 213
EXPO = CPBYCV / Y1 TO = T1 * (1.0 + Y2) RETURN END SUBROUTINE CONTRL (Ml , PRl , TEMP , AREA, PR2 ) C THIS SUBROUTINE CALCULATES THE INLET STAGNATION PRESSURE NEEDED C FOR OBTAINING A GIVEN MASS INFLOW. REAL*8 M1,MACH,MSQ, S1,S2,S3,S4,S5,PR1,PR2,RUNIV,RFUEL,MEANCP, 1 CPBE, CPBYCV, SBCONS,EMISS, AREA, TEMP COMMON /CONST/ RUNIV, RFUEL, CPBYCV, MEANCP , CPBE , SBCONS , EMISS 51 = DSQRT (RFUEL * TEMP / CPBYCV) 52 = Ml / (AREA * PRl * 1.013E+5) MACH = SI * S2 MSQ = MACH* *2 53 = CPBYCV 1.0 54 = CPBYCV / S3 55 = 1.0 + (S3 * MSQ / 2.0) PR2 = PRl * (S5**S4) RETURN END SUBROUTINE INTEMP (TEX, T04, TOl, TRAD, T02, DELTAT) C THIS SUBROUTINE CALCULATES THE CORE INLET STAGNATION TEMPERATURE C GIVEN THE CORE EXIT TEMPERATURE AND BALANCE OF PLANT CONSTANTS C LIKE TURBINE EFFICIENCY, COMPRESSOR EFFICIENCY, TURBINE ENTHALPY C EXTRACTION RATIO ETC. . INTEGER N REAL*8 TEMPS (2) , TEX, TOl, T02, T04, ETAT, ETAC, EPSIL, UBARA, TERM, 1 DELTAT ETAT =0.75 EPSIL =0.20
PAGE 214
206 ETAC =0.8 UBARA = 1.8E+5 TERM = 1 (ETAT / (ETAT EPSIL) ) TEMPS (1) = TOl TEMPS (2) = TRAD T04 = TEX * (1 EPSIL) CALL RK42 (TEMPS, T04) TOl = TEMPS (1) TRAD = TEMPS (2) T02 = TOl * (1 TERM / ETAC) DELTAT = T02 TOl RETURN END SUBROUTINE TEMP (TIN, MIN, MEX, STEP, AREAl , BEMASS, POWER, MASS, 1 VEL,COREPR,TEX,CORTMP,BETMP,HEATNG) C THIS SUBROUTINE CALCULATES THE CORE FUEL TEMPERATURE AND EXTERNAL C BERYLLIUM TEMPERATURE BASED ON POWER GENERATED WITHIN THE CORE C PLUS RADIATIVE AND CONVECTIVE HEAT TRANSFER FROM THE CORE TO THE C MODERATOR. IT ALSO ACCOUNTS FOR THE GAMMA AND NEUTRON HEATING OF C THE MODERATOR. REAL*8 RENOLD, PRANDL, RUNIV, RFUEL, CPBYCV, MEANCP, CPBE, SBCONS, 1 EMISS, VOL, HEIGHT, DIA, RADIUS, AREAIN,AREAEX, THROAT, 1 BETHIK, ENRICH, TIN, STEP, AREAl , BEMASS, POWER, MASS, MACHSQ, 1 VEL , COREPR, CORTMP , BETMP , HEATNG ( 3 ) , TERM ( 7 ) , MOLWT , 1 RHO,MU,CONDUT,MIN,MEX,TEX COMMON /CONST/ RUNIV, RFUEL, CPBYCV, MEANCP, CPBE, SBCONS, EMISS COMMON /CORDIA/ VOL, HEIGHT, DIA, RADIUS, AREAIN,AREAEX, THROAT, 1 BETHIK, ENRICH
PAGE 215
207 PRANDL = 0.2355 MOLWT =28.2 RHO = COREPR * MOLWT / (0.08205 * CORTMP) MU = 3.640E8 * CORTMP + 2.27E5 RENOLD = RHO * VEL * DIA / MU CONDUT = 1.810E4 * CORTMP + 0.102 TERM(l) = POWER * 0.9 HEATNG(l) = POWER * 0.1 TERM(2) = AREAl * SBCONS * EMISS * (CORTMP * 4 BETMP ** 4) HEATNG(2) = TERM (2) TERM(3) = CONDUT * 0.0364 * (RENOLD * PRANDL ) ** 0.75 / DIA TERM(5) = AREAl * TERM(3) * (CORTMP BETMP) HEATNG(3) = TERM(5) TERM(6) = HEATNG(l) + HEATNG(2) + HEATNG(3) TERM(7) = MEANCP * ( (MEX * TEX) (MIN * TIN)) TEX = TEX + (STEP / (MASS * MEANCP)) * (TERM(l) TERM(2) 1 TERM(5) TERM(7)) CORTMP = (TEX + TIN) / 2 C BETMP = BETMP + (STEP / (BEMASS * CPBE) ) * (TERM(6) * 0.7) BETMP = BETMP + (STEP / (BEMASS * CPBE)) * (TERM(6) * 1.0) RETURN END BLOCK DATA REAL*8 TOTALl (5) ,TOTAL2 (5) ,FISS1 (5) ,FISS2 (5) ,ABS1 (5) ,ABS2 (5) , 1 NUl (2) ,NU2 (2) ,TRANS1 (5) ,TRANS2 (5) ,ONET01 (5) ,ONET02 (5) , 1 TWOT02 (5) ,D1 (2) ,D2 (2) ,SIGA1 (2) ,SIGA2 (2) , SIGFl (2) , 1 SIGF2 (2) , SIGH (2) , SIG12 (2) , SIG22 (2) , SIGTRl (2) , SIGTR2, SIGRl (2),SIGR2(2) ,NU(2),RHO(14) ,MASS(14) ,H(7) ,SECOND(8) , LSTAR,BETA(6) ,ALAMDA ( 6) , VOL, AREAIN,AREAEX, THROAT, RADIUS, HEIGHT, DIA, BETHIK, ENRICH, RUNIV, RFUEL, CPBYCV, MEANCP, TIME, TIMAX, TLOOP, STEP, SBCONS, CPBE, EMISS 1 1 1 1 1 COMMON /MICRO/ TOTALl, TOTAL2, FISSl, FISS2,ABS1,ABS2, NUl, NU2, TRANSl , TRANS2 , ONETOl , ONET02 , TWOT02
PAGE 216
208 COMMON /MACRO/ Dl, D2, SIGAl, SIGA2, SIGFl, SIGF2, SIGH, SIG12, SIG22, 1 SIGTR1,SIGTR2,SIGR1,SIGR2,NU COMMON /INTPRM / LSTAR, BETA, ALAMDA COMMON /CORDIA/ VOL, HEIGHT, DIA, RADIUS, AREAIN, AREAEX, THROAT, 1 BETHIK, ENRICH CCMMON /CONST/ RUNIV, RFUEL, CPBYCV,MEANCP, CPBE, SBCONS, EMISS COMMON /TIMECO/ TIME, TIMAX, TLOOP, STEP DATA TOTALl,TOTAL2 /21 . 64, 23 . 39, 4 . 50, 1 . 804, 5 . 29, 1 291.52,10.16,4.11,0.845,6.199/ DATA FISS1,FISS2 /8 . 02, 8 . 27E2, 3*0 . 0, 233 . 58, 4*0 . 0/ DATA ABS1,ABS2 /12 . 31, 7 . 37, 5 . 02E3, 5 . 682E5, 4 . 125E3, 1 276.8,1.244,4.265E3,3.107E3,4.894E3/ DATA TRANS1,TRANS2 /7 . 901, 8 . 175, 3 . 173, 2 . 633, 3 . 716, 1 291.52,10.16,3.971,0.2844,6.015/ DATA ONETOl /9 . 324, 16 . 02, 4 . 48, 1 . 794, 5 . 219/ DATA ONET02 /3 . 36E3, 2 . 376E3, 1 . 331E2, 1 , OlE2, 7 . 89E2/ DATA TWOT02 /14 . 717, 8 . 917, 4 . 108, 0 . 842, 6 . 194/ DATA NU1,NU2 /2 . 431, 2 . 781, 2 . 419, 2 . 320/ DATA LSTAR /6.0E4/ DATA BETA /O . 000215, 0 . 00142, 0 . 00127, 0 . 00257, 0 . 00075, 0 . 00027/ DATA ALAMDA /O . 0124, 0 . 0305, 0 . Ill, 0 . 301 , 1 . 1, 3 . 0/ DATA VOL, HEIGHT, DIA, RADIUS, AREAIN, AREAEX, THROAT, BETHIK, 1 ENRICH/1.57,2.0,1.0,72.1,0.09,0.09,0.05,50.0,0.85/ DATA RUNIV, RFUEL, CPBYCV,MEANCP, CPBE, SBCONS, EMISS /8.314, 1 294. 98, 1.403, 1000. 0,2428. 46, 5. 67E8, 0.25/ DATA TIME, TIMAX, TLOOP, STEP /O . 0, 10 . 0, 0 . 1, 0 . 0001/ END
PAGE 217
APPENDIX B BPGCR SYSTEM EQUATIONS. AND THE DERIVATION OF THE BPGCR TRANSFER FUNCTIONS Chapter 3 of this dissertation describes the transition analysis of the BPGCR from a low power mode of operation to a high power mode of operation. In Chapter 3, brief descriptions of the equations used for modelling the system are given. This Appendix describes the full set of equations as used in Chapter 3 and Chapter 4 of this dissertation. Chapter 4 includes a linear stability analysis of the BPGCR. For this analysis, the nonlinear system of equations are linearized and transfer functions are derived from this set of linearized equations. This Appendix describes the procedures used for linearizing the system equations, and gives the methodology for deriving the BPGCR transfer functions. It also contains the transfer functions used in Chapter 4. A list of variables used in this Appendix is given at the end of this Appendix. 209
PAGE 218
Relationship Between Core Inlet Temperature. Averag e Core Temperature. Core Power and Core Exit Temperature 210 The equations relating the average core fuel temperature and the average beryllium moderator temperature are given by Mf (t) = N(t) h Aw (Tf Tbc) C^ (m3'(t) Ts' m 2 (t) T 2 ) (B1) and. Mbc C?" = h Aw (Tf Tee) Ci where " Ci" represents a constant heat removal term from the (B2) beryllium moderator, i.e., it is assumed that a constant amount of heat is removed from the beryllium moderator/ reflector at all times, irrespective of the core power or other operating conditions. This could be accomplished, for example, by passing a gas such as helium through the beryllium wall at variable flow rates. All other variables mentioned in the above equations are described at the end of this Appendix. The state points 2' and 3' used in the Equation (B1) represents the conditions inside the BPGCR core just after the inlet nozzle and just before the exit nozzle, respectively. Figure B1 shows the state points (1, 2, 3, etc.) used elsewhere in this Appendix, but it does not show the points 2' and 3'.
PAGE 219
211 CO Figure B1 . TS Diagram for tÂ±ie BPGCR Power System
PAGE 220
212 And finally, the core inlet temperature, the core exit temperature and the core (fuel) average temperature are related by the expression Tf = T2' + T3' 2 (B3) For performing a linear stability analysis of the BPGCR, we have to derive the transfer functions relating the variables involved. To obtain the transfer functions, the above set of equations have to be linearized. The procedure used is as follows: first all the variables in the above equations are written as a sum of two partsÂ— one part being their steady state value and the other being their variations from their steady state values. As an example, the core thermal power, N(t), is rewritten as N(t) = Nss + 3N(t) where 'Nss' denotes the steady state part and '9N' denotes the time dependent variation from the steady state value (in a linear stability analysis, it will assumed that the magnitude of the time dependent part is small compared to the steady state value). All other variables are similarly written as a sum of a steady state value and a time dependent variation from the steady state value. Equations (B1), (B2), and (B3) are rewritten in terms of these new variables. Subtracting the steady state equations from this new set of equations results in a third set of equations given belowÂ— which involves only the
PAGE 221
213 derivatives of the time dependent part of all the variables given in Equations (B1). (B2), and (B3). Mp = 3N(t) h Aw OTf 3TBe) [mÂ«(3T3' BT 2 ) + Ti am3' 1*2^ am2] (B 4) Mse C?" = h Aw OTf 9 TBe) (5.5) aTf = 3T2 + 9T3' 2 (B6) In arriving at Eq. (B4), the following assumptions have been made; a) The fuel mass inside the core is taken to be a constant, M(ss, for estimating the fuel heat capacity term (i.e., MfCpf term on the right hand side of Eq, (B4)). After the system becomes critical, the fuel mass inside the core varies very little for any significant reactivity variations (see, e.g.. Fig, 3.19). This is due to the fact that the core reactivity is a strong function of fuel mass inside the core and, hence, the total fuel mass inside the core varies only slightly from its steady state value for the kind of reactivity variations allowed by a linear stability analysis. (Note, however, that it has not been assumed that core inlet and outlet fuel mass flow rates are the same or they remain constant to provide a constant fuel mass inside the core.) b) Terms involving the product of two small variations from the steady state values (i.e., the "9" quantiies) are neglected in comparison to other terms.
PAGE 222
214 Now to derive the transfer functions from the above set of equations, the Laplace transform of Eqs. (B4), (B5) and (B6) are taken to obtain Eqs. (B7), (B8) and (B9) given below; _ P s 9 Tf(s) = =^i^N(s) 9 Tf (s) + aiBe(s) h h A\v ^[mOT3 (s) 3 T 2 (s)) + Tf 8013 (5) 8012(8)] hAw (B7) Â— s8TBe(s) = 8ff(s) 8TBe(s) (B8) hAw afrfs) = 2 * (B9) These equations can be rewritten as ai s 3Tf{s) = a 2 3N(s) 9Tf (s) + 9TBe(s) as 3T2(s)] as Ts^ 8oisÂ’ + as l2^ 8012 (B10) a4 s 8Tbc(s) = 8Tf(s) 8Tbc(s) (B11) afKs) = 2 (B12) where mP ^ = = ^ h A^ , a2 = Â— h A^ , as = h A^ , and and
PAGE 223
Equation (B11 ) can be solved for 8 TBe to obtain ^Be ^ 1 3Tf a 4 S+l . (B13) (The right hand side of Eq. (B13) is denoted as the transfer function G 8 in Chapter 4. Similarly. Eq. (B12) provides an expression for the transfer function G7 used in Chapter 4, viz., G7 = 1/2.) Using Eq. (B13) to eliminate 5TBe from Equation (B10) and using the relationshiop (B12), we can finally obtain the relationship between the core inlet temperature, the core exit temperature, the core inlet and exit fuel mass flow rates and the core power as as = a 68 T 2 + a20N + a3TÂ®3m2 asTf^dma (B14) where ^ ai a4S _ ^5 = T ^ ^ 2(1 ^ a 4 s) ^ ^3Â”^ , and aÂ« = ^s^ ^ 2 2(1 + a 4 s) asmf Equation (B14) can be rewritten in a form compatible with the transfer function notation used in Chapter 4 as ^ 3 ' ~ Gi 0T2' + G23N + G39m2' + G43m2' (B15) where
PAGE 224
216 G2 = ^ ^ as G3 = ^ /X'SS as T2 as , and G4 = as T2' as Balanceof Plant Model Equations In the BPGCR system being analyzed in this dissertation, there is a converging nozzle at the core inlet and a convergingdiverging nozzle at the core exit (see Fig. 33). It is assumed in this analysis that the flow through these nozzles is isentropic and, hence, there is no stagnation temperature (or pressure) change in going through these nozzles. The temperatures T 2 ' and TaÂ’ used in Equations (B1) through (B15) represents the conditions inside the core Just after the inlet nozzle and just before the exit nozzle, respectively. The points 2 and 3 shown in Fig. B1 represents conditions outside the core, just before the inlet nozzle and just after the exit nozzle, respectively. Because of the isentropic flow assumption, we can write T 2 Â’ = T 02 ' = To 2 , and 3 = i03 = i03 where the subscript "0" refers to stagnation conditions. When the fuel velocity is small (a few tenths of a Mach number), the total temperature and the stagnation temperature can be interchanged with very little error.
PAGE 225
217 Turbine Equations In the actual BPGCR system, an MHD generator is employed to convert the thermal power generated Avithin the BPGCR to electrical power. For the purpose of this dissertation, this MHD generator is modelled as a conventional turbine (see Chapter 3 for more discussion on this). So, the nozzle and diffuserÂ— normally associated with an MHD generator, are not included here. For a turbine, the pressure ratio, rt is given by rt = Po3/Po4 (B16) (see Fig. B1 for the location of the state points 3 and 4). The turbine enthalpy extraction ratio, e, and isentropic efficiency, rit, are defined as e = (To3 To4) / To 3 . and (B1 7) Â•It = (Tq3 T04) / (Tq3 To4i) . (B1 8) Using Equation (B1 7) and the isentropic relationship between pressure and temperature, viz.. we get Tq3 _ Tq41 (B19) To4 = Tq 3 (1e), and
PAGE 226
218 JIlll To4t = [rt]\ Y I To3 . Using the above two relations, we can obtain the expression for the isentropic turbine efficiency as ^t = [1m Y and, hence. rt = Tit Tite niJL Yl So, given the turbine inlet conditions and the turbine parameters like isentropic efficiency, enthalpy extraction ratio, etc., the conditions at the turbine exit can be obtained from the relations Tq4 = Tq 3 (1 e) (B20) and Po4 = Po3 Tit e (B21) The equation for the turbine exit temperature, Tq 4 , (i.e., Eq. B20) can be linearized as follows 3To 4 = (1e) 3 Tq3 Â• Or, using the Laplace transform variable. dTo4 (s) = (1e) dTTos (s) . (B22)
PAGE 227
219 Heat Exchanger/Radiator Equations After exiting from the turbine, the fuel/working fluid (gas mixture) passes through a heat exchanger to transfer the waste heat to a radiator. It is assumed, for this analysis, that an average heat exchanger temperature and an average radiator temperature can be used to describe the heat transfer between the two. With this assumption, the relevent equations for the average heat exchanger and radiator temperatures are as given below; [Mf^ cfp + ^ cfp(T04 To ,) t^hx ^hx ("I^hx ^rad) (B23) [Msea, ^etal ^ = ^hx ^hx Trad) " ^rad Â£ Trad and. (B24) Thx = Tqi + Tq 4 2 (B25) As was done for the core fuel temperature equations (Equations (B1), (B2) and (B3) ), the Equations (B23), (B24) and (B25) can be first linearized by considering small variations around some mean (i.e., rewriting all the variables as a sum of a steady state term and a time dependent variation around this steady state value) and then subtracting the steady state equation from this new equations. Finally, the Laplace transform of the equations are taken to obtain the following set of equations
PAGE 228
220 [S] 5 Th, = (5to4 8To i) (sThx 5 T,ad) STrad (S) (C 4 +S)= C 3 5fhx (S) 5 Thx (s) = + 5 Tq4(s) 2 where (B26) (B27) (B28) C 3 = ^hx ^hx M^IFai + Mhe . and C 4 = C 2 + C 5 p ^Â»^ne ''p aT^ad' .aTradJss with C 5 = ^ ^ M^fFai + Mhe . (Note that before Laplace transform of Eq. (B24) can be taken, it has to be linearized because of the nonlinear term involving Trad"^) The three Equations, (B26), (B27), and (B28) can be solved simultaneously to finally obtain an expression for the temperature at the outllet of the heat exchanger (Toi) in terms of the heat exchanger inlet temperature (T04) as follows 5Toi = 2 Ce C7 + 2 C0 +C7 C2 C7 (C 4 + s) C2C7 (C 4 + s) 5To4 (B29) where
PAGE 229
221 C2 . and Q_ _ hhx Ahx C2 or, taking the Laplace transform of Eq. (B29), we can write it as Compressor Equations The isentropic efficiency of the compressor, ric. is defined as "He = (Toi T021) / (Toi T02} (see Fig. B1 for the location of the state points 1 and 2) so. The isentropic relationship between pressure and temperature is given by The pressure ratio on the right hand side of Eq. (B31) is known, once the turbine parameters (efficiency, enthalpy extraction ratio, etc.) are given. From the last two relations given above, we get 6 Toi(s) = Cs 5 To 4 (s) . (B30) r lYl To2i _ Pq2 ~T~ Toi I.P01. (B31) ^02 ^ 1 . Y1 TOI ^c Y
PAGE 230
222 or, in terms of Laplace transform variable. dTp2(s) dToi(s) 1^ 1Tl Y ^c L iPoii (B32) Equations (B22), (B30) and (B32) can be combined to obtain the complete transfer function of the balance ofplant, except for the time delay functions, G5, in Chapter 4 as STq2 (s) 5To3(s) = (Cs) (1 e) Tite) (B33) Time Delay As mentioned in Chapter 4, the total delay of the external loop is divided into two parts; one for the transport of fuel from the core to the heat exchanger and the other for the transport of fuel gas mixture from the heat exchanger back to the core. We can write an equation for the time delay of the form To4(t) = To3(t ti ) where 'xi' time delay for the transport of the fuel from the core to the heat exchanger and T 03 ' and T 04 ' are the core exit temperature and the heat exchanger inlet temperature, respectively. Taking the Laplace transform of the above equation gives 5 To 4 (s) = 3 To 3 exp(xis) which can be approximated by (for small values of xi) as
PAGE 231
223 3To4 (s) = aTo3(s) / (1 +Ti s) . (B34) (The typical values of the external loop circulation times used in this dissertation are of the order of a tenth of a second. Since, "xi " is only half of this, the above approximation should not be too bad.) Similarly, 3To2 (s) = 3Toi (s) / (1 + X 2 s) where "X 2 " is the time delay for transporting the fuel from the heat exchanger back to the core. Equations for the Core Pressure, the Core Inlet and Exit Fuel Mass Flow Rates and the Core fuel Mass Core Average Pressure The core average pressure is given by the relation Pdt) = n(t) TKt) Ru e / V where "e" is a constant for conversion of pressure from atmosphers to pascals (see the end of the Appendix for a list of symbols for the other variables). This equation can be written in terms of the fuel mixture (UF6 + He) mass inside the core as Pf (t) = ei Mf(t) Tf(t) (B35) where "ei" is a constant. For small deviations for steady state, Eq. (B35) can be linearized as follows
PAGE 232
5P(t) = ej [ffÂ® 5M(t) + MfÂ® 8 Tf (t)] or, taking the Laplace transform of the above equation, 3P(s) = C 2 3Mf{s) + C 3 9Tf(s) 224 (B36) where "e 2 " and "ea" are constants. Core Fuel Mass The change in the core fuel mass occurs due to a difference in the core inlet and exit fuel mass flow rates, i.e., (dMf/dt) = m 2 '(t) ms'lt) (under steady state conditions, m 2 = 013 '= mfss ) Considering small variations from steady state, we can write the above equation as (d3M/ dt) = 9 m 2 ' 9 m 3 ' or, taking the Laplace transforms of the above equation, s 9M({s) = 9 m 2 ' (s) 9m3' (s) or, 9Mf = (1/s) (9m2Â’ 9ni3') . (B37) (In Chapter 4. the transfer functions G13 and G14 are given by G13 =l/s and G14 = 1/s.) Core Inlet Fuel Mass Flow Rate The fuel mass flow through the inlet nozzle is given by ni 2 '(t) = d 2 '(t) V 2 ' (t) Am (B38)
PAGE 233
225 with d2' given by d2'(t) = P(t) / (RfT2Â’(t)) . In the above equation, as an approximation, the core inlet pressure is taken to be the same as the core average fuel temperature, and V2' (t) = (yRfT2')i/2 M2'(t) . The core inlet Mach number, M2', can be written in terms of the core average pressure and inlet stagnation temperature, and after simplifying the above eqation, (B38), we finally get am2'(t) = fi 3P(t) + f2 aTo2 (t) where and Cg is a constant. In terms of the Laplace transform, we can write am2'(s) = fi ap(s) + f2 aTo2 (s) . (B39) Core Exit Fuel Mass Flow Rate At the core exit, the fuel mass flow is assumed to always be choked. So, we can write the expression for the maximum flow through a convergingdiverging nozzle (which occurs when the nozzle is choked) as
PAGE 234
226 riig(t) = Pqc Aex VRfToc Â• (B40) The core stagnation pressure is assumed to be the same as the core average pressure and the core stagnation temperture is taken to be equal to the the core average temperature. Since the fuel velocity within the core is very small (only around 30 m/s or so), these assumptions should be fairly good. Linearizing the above equation, we get am3'(t) = gi aXf+ g2 3P where gi = C Pss (1/2) (Tss)3/2 . g 2 = C (Tss)(i/2) , and "C" is a constant. Taking the Laplace transform of the above equation results in the following equation for the core exit fuel mass flow rate 3m3'(s) gi aTf(s) + g 2 aP(s) . (B41)
PAGE 235
227 List of Variables Used The following is a list of variables used in this Appendix: Aex = Throat area of the convergingdiverging nozzle at the core exit (m2), Ahx = Effective heat transfer area of the heat exchanger primary (m2) Arad = Surface area of the radiator (m2). Aw = Wall area of the core (m2) , CpBe = Average specific heat at constant pressure of the beryllium used as the reflector /moderator (J/KgK), Cpf = Average specific heat at constant pressure of the fuel (UFe+He) mixture (J/KgK). Cphe = Average specific heat at constant pressure of the helium used in the secondary side of the heat exchanger (J/KgK), CpHietai = Average specific heat at constant pressure of the metal used in the heat exchanger primary/ secondary side (J/KgK), h = Average convective heat transfer coefficient for the fuel/beiyllium wall heat exchange inside the core (w/m2K), hhx = Average convective heat transfer coefficient for the heat transfer between the primary and the secondary side of the heat exchanger (w/m2K), M = Mach number of the fuel,
PAGE 236
228 Mbc = Mass of the beryllium used as the moderator/reflector (Kg), MfC = Fuel mixture (UFe+He) mass within the core (Kg), MfÂ«s = Steady state fuel mass inside the core (Kg), Mfiix = Mass of the fuel within the primary side of the heat exchanger (Kg), Mhe = Mass of the helium within the secondary side of the heat exchanger (Kg), MmetaiP^ = Mass of the metal used (e.g., as tubes) in the primary side of the heat exchanger (Kg), MmetaiÂ®Â®^ = Mass of the metal used (e.g., as tubes) in the secondary side of the heat exchanger (Kg), m2= Fuel mass flow rate at the core inlet (Kg/sec), ms= Fuel mass flow rate at the core exit (Kg/sec), mss = Steady state fuel mass flow rate into (or from) the core (Kg/sec), mftix = Average fuel mass flow rate through the heat exchanger (Kg/sec), N = Core thermal power (watts), n = Number of moles of fuel mixture (UFe+He) within the core, P = Pressure variable (Pascals), T2= Core inlet fuel temperature (K),
PAGE 237
229 Ta= Core exit fuel temperature (K), Tf = Core average fuel temperature (K), Tbc = Average beryllium temperature (K), Toi = Fuel temperature at the heat exchanger outlet (K), Tq 4 = Fuel temperature at the heat exchanger inlet (K), Thx = Average fuel temperature in the primary side of the heat exchanger (K), Trad = Average temperature of the radiator (K), t = time variable, Rf= Fuel mixture gas constant (Nm/KgK), s = Laplace transform variable, V = Core volume (m3) , V = velocity of the fuel (m/sec), V = Ratio of specific heats for the fuel mixture, O = Stefan Boltzmann constant (w/m2K4), and Â£ = Emissivity of the radiator material Note: A " " on the top of variables indicate average or mean quantities.
PAGE 238
LIST OF REFERENCES 1. G.I. Bell, "Calculation of the Critical Mass of UF6 as a Gaseous Core, with Reflectors of D 2 O, Be aind C," USAEC Report LA1874, Los Alamos National Laboratory, Los Alamos, New Mexico (February, 1955) 2. J,A. Angelo, Jr., and David Buden, "Space Nuclear Power," Orbit Book Company, Malabar, Florida (1985) 3. M.M. Panicker, "Neutronics of a Coupled Multiple Chamber Gaseous Core Reactor Power System." Ph.D Dissertation. University of Florida, Gainesville, Florida (April, 1989) 4. F.E. Rom. "Gaseous Nuclear Rocket," Patent No, 3,202,582 (August, 1961) 5. R.G, Ragsdale, "Are Gas Core Nuclear Rockets Attainable," NASA TM X524 (1968) 6. R.T, Schneider, and K. Thom, Editors, "Proceedings of a Symposium on Research on Uranium Plasmas and their Technological Applications," Gainesville, Florida (January, 19701 7. T.S. Lantham, "Nuclear Studies of the Nuclear Light Bulb Rocket Engine," UARL Report No. G9 103753, United Aircraft Research Laboratories, East Hartford, Connecticut (September, 1968) 8. E.C, Gritton, and B. Pinkel, 'The Feasibility of the Gaseous Core Nuclear Reactor Concept for Electric Power Generation," Rand Corporation Report RM5721PR, The Rand Corporation, Santa Monica, California (June, 1969) 9. J.R, Williams, and S.V. Shelton, "Gas Core Reactors for MHD Power Systems," "Proceedings of the Symposium on the Research on Uranium Plasmas and their Technological Applications," Gainesville, Florida (January, 1970) 230
PAGE 239
10 . R.T. Schneider, and M.J. Ohanian, "Nuclear Engine," Patent Disclosure to NASA (January, 1970) 231 11. N. J. Diaz, E.T. Dugan, and C.C. Oliver, "Neutronics and Energetics of Pulsed Gaseous Core Nuclear Systems," Final Report to NSF Grant No. ENG 7501437, University of Florida, Gainesville, Florida (April, 1978) 12. E.T. Dugan, "The Nuclear Piston Engine and Pulsed Gaseous Core Power Systems," Ph.D Dissertation, University of Florida, Gainesville, Florida (March, 1976) 13. N.J. Diaz, and E.T. Dugan, "Gas Fueled Heterogeneous Core Reactors," Disclosure of Invention, University of Florida, Gainesville, Florida (October, 1974) 14. K.l. Han, "Heterogeneous Gas Core Reactor," Ph.D Dissertation, University of Florida, Gainesville, Florida (December, 1977) 15. INSPl Quarterly Progress Report, lNSPlQRUF86001, U.S. Air Force Contract No. DNA00185C0329, University of Florida, Gainesville, Florida (March, 1986) 16. P.M, Gresho, J.E. Bickel, L.R. James, and J.S. McDonald, "SNAP2 Startup Analysis," Atomics International Report No. NAASR9078, Atomics International, Canoga Park, California (January, 1964) 17. S. Birken, "Response of SNAP8 Reactor During Automatic Startup," Atomics International Report No, NAASR9646, Atomics International, Canoga Park, California (Sept., 1964) 18. T.S. Lantham, H.E. Bauer, and R.J. Rodgers, "Studies of Nuclear Light Bulb Startup Conditions and Engine D)mamics," UARL Report No. H9 103754, United Aircraft Research Laboratories, East Hartford, Connecticut (September, 1969) 19. T.S. Lantham, H.E. Bauer, and R.J. Rodgers, "Analytical Studies of Startup and Dynamic Response Characteristics of the Nuclear Light Bulb Engine," UARL Report No. J9109005, United Aircraft Research Laboratories, East Hartford, Connecticut (September, 1970) 20. L.M. Petrie, and N.M. Greene, "XSDRNPM: AMPX Module with One Dimensional Sn Capability for Spatial Weighting," AMPX: A Modular Code System for Generating Coupled Multigroup
PAGE 240
232 NeutronGamma Libraries from ENDF/B, Oak Ridge National Laboratory, Oak Ridge, Tennessee (1976) 21. J.F. Breismeister, Ed,, "MCNPA General Monte Carlo Code for Neutron and Photon Transport, Version 3A," LA7396, Los Alamos National Laboratory, Los Alamos, New Mexico (1986) 22. F.P, Incropera, and D.P. Dewitt, "Fundamentals of Heat Transfer," John Wiley and Sons, New York (1981) 23. A. H. Shapiro, "The Dynamics and Thermodynamics of Compressible Fluid Flow," Vol, 1, Ronald Press Co., New York (1953) 24. 1, K, Kikoin, V. A. Dmitrievsky, Y. Y. Glazkov, 1. S. Grigoriev, B. G. Bubovsty, and S. V, Kersnovsky, "Experimental Reactor with Gaseous Fissionable Substance (UF6)," "Proceedings of the Second International Conference on the Peaceful Uses of Atomic Energy," Geneva, Switzerland, Vol. 2, p. 232 (1959) 25. E. T. Dugan, N. J. Diaz, and E. E. Carroll, Jr., "Gas Core Reactor Neutronics Theoretical Modeling and Experimental Verification," Nuclear Technology, Vol, 69, pp. 134153 (1985) 26. D. L. Hetrick, "Dynamics of Nuclear Reactors," The University of Chicago Press, Chicago (1971) 27. J. C. Vigil, "ANCON User's Manual," LA4616, VC32, TlD4500, Los Alamos National Laboratory, Los Alamos, New Mexico (1958) 28. G.R. Keepin, "Physics of Nuclear Kinetics," Addison Wesley Publishing Company, Reading, Massachusetts (1965) 29. M. A. Schultz, "Control of Nuclear Reactors and Power Plants," McGrawHill Book Company, New York (1961) 30. J. MacPhee, 'The Kinetics of Circulating Fuel Reactors," Nuclear Science and Engineering, Vol, 4, No. 4, p. 588 (1958)
PAGE 241
BIOGRAPHICAL SKETCH Kiratadas Kutikkad was bom on August 26, 1958, in Kerala, India. After graduating from high school in 1974, he joined the Government Victoria College affiliated with the University of Calicut. He obtained a Bachelor of Science degree in physics in 1979 and a Master of Science degree in applied physics in 1981. The same year, he entered the graduate program at the Indian Institute of Technology, Kanpur, India, and graduated from there with a Master of Technology degree in nuclear engineering in 1984. After working in the same institute as a research engineer for six months, he obtained admission from the University of Florida, Gainesville, and started working towards his Ph.D, in nuclear engineering sciences. As a graduate student at the University of Florida, he has worked as a teaching assistant in the Nuclear Engineering Department as well as in the Physics Department. He has also worked as a graduate research assistant for the Innovative Nuclear Space Power Institute at the University of Florida. He has been working as a research scientist at the University of Missouri, Research Reactor Facility, since January, 1990. 233
PAGE 242
I certify that I have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. Edward T. Dugan, Chairman Associate Professor of Nuclear Engineering Sciences I certify that 1 have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. Alan M. Jacc^, Professor of Nuclear Engineering Sciences 1 certify that 1 have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. Nils ~P Professor of Nuclear Engineering Sciences
PAGE 243
I certify that I have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. Calvin C. Oliver Professor of Mechanical Engineering 1 certify that 1 have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. tobert J. H: Professor o L rahan Chemistry This dissertation was submitted to the Graduate Faculty of the College of Engineering and to the Graduate School and was accepted as partial fulfillment of the requirements for the degree of Doctor of Philosophy. December, 1991 M. Phimps 'ollege 01 Engineering Madelyn M. Lockhart Dean, Graduate School
