ML20115G551
| ML20115G551 | |
| Person / Time | |
|---|---|
| Site: | 05000128 |
| Issue date: | 08/31/1995 |
| From: | Parish T, Rearden B TEXAS A&M UNIV., COLLEGE STATION, TX |
| To: | |
| Shared Package | |
| ML20115G536 | List: |
| References | |
| NUDOCS 9607190198 | |
| Download: ML20115G551 (138) | |
Text
{{#Wiki_filter:iii ABSTRACT Engineering Analysis of a Power Upgrade for the Texas A&M Nuclear Science Center Reactor. (August 1995)
'Bradley Thomas Rearden, B.S., Texas A&M University Chair of Advisory Committee: Dr. Theodore Parish Heat transfer, neutronics, and accident analyses are presented to support an increase of the power level of the Texas A&M University Nuclear Science Center Reactor. The upgraded steady state power level is to be 1.49 MW, from the current level of 1.0 MW.
Several computer codes were used to perform the safety analyses needed to support a licensing change. These codes were WIMSD4/m, DIF3D, NCTRIGA, and PARET. The I limiting reactor parameters were calculated, and standard reactor accident scenarios were analyzed. The results of the analyses demonstrated that the power upgrade could be achieved without creating any significant safety hazards or restricting operation beyond the current limits. Included in this thesis is a description of the relevant reactor systems and verification of the numerical methods. l l 1 1 1 1 . l l 9607190198 960715 PDR ADOCK 05000128 P PDR
i iv h ACKNOWLEDGMENTS I would like to.thank the members of my committee, Dr. Lazarov, Dr. Parish, and ; Dr. Reece, for all of their help in completing this project. Dr. Parish deserves a special l i thanks for his help in getting the codes to run when it seemed impossible and for checking l my work, time after time, to make sure my presentation was as good a possible. Thanks to ! the Nuclear Science Center for supporting me over the past months while working on this i project. The operations staff provided valuable help in understanding the systems and gathering experimental data. I would also like to thank Chien-Hsiang Chen for his help in ! drawing several of the figures included in this thesis. l 1 Thanks to my parents, Tom and Sue, for their love and encouragement through my l many years here at A&M. And, thanks to Kristin for her understanding while I've worked long hours to complete this project. l l 1 i 1 l l i
1 1 I ! y I l TABLE OF CONTENTS Page ABSTRACT. .. . . . .. .. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .. . ...iii i A C KN OWL ED G M ENTS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .. . . . .. . .. . . . . . . . . .. . . . . . . . . . . . . . iv l TAB L E O F C ONTENTS . . . . . . . .. . . . . . . .. .. . . .. . . . .. . . ... . . . . .. . . . . . . . . . . . .. . . .. . . . . . . . ........v l l LIST OF FIGURES. . . . . . ..........................................................vii i l LIST OF TABLES.. . ..................................................................xi l l CHAPTER I INTRODUCTION. . . . . . . . . . ................................................I 11 REACTOR DESIGN.. . .. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ......5 Mechanical Design . . . . .. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .........5 Nuclear Design .. .. . . . . . . . . . . . . . . . . . . . . . . . . . . . .. . . . . 15 Fuel Description and Safety Limits. .... ...... . .. ... .. . . . . . . . . 18 i l Thermal Design . .. . .. . . . . . . . . . . . . . . . . . . . . . . . . . . ... . 22 t III NUMERICAL ANALYSIS TECHNIQUES . .... ... ..... ........ .. . . . . . . . . . . . . . . . . 24 Neutronic Analysis .......... ....... . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .. 24 i i Thermal Hydraulic Analysis .. . . ... ...... ......... . . . ..... . ........ .. . ... .... 3 2 < l , i Core Average Temperature... ...................................................49 1 i IV S AFETY EVALU ATIO N .. . ... ... . . . . . .. .... ................ .... .. .... .. .. ...... . . . 50 R eactor Analysis . . . . . . . .. . .. . . . . . . .. . . . . .. . .. . . . . .. . . . . . ... . .. . . . . . . . . . . . . . .. . 5 0 Evaluation of the Limiting Safety System Setting ............ .. . .. ... ... .... 59 LSSS for Steady State Operation... . .. . . . . . . . . . ................60 The LSSS and Maximum Allowable Reactivity Insertion for Pulsing.. . 63 l Accident Analysis..... . .. .. ...................................................70 1 V CONCLUSIONS.... . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 83 Future Work .. ... .. ........ .... . . . . . . . . . . . . . . . . . . . . .. . . . . 8 5 REFERENCES. . .. . . . . .................................................86 APPENDIX A - WIMSD4/M MODELING.. . .... . . . . . . . . . ......................88 Fuel Rod Model ... ............... .. ..........................................88 Control Rod Model .. ....... .. . ...... ....... . . . .......................89
I > t
,b vi i k
l, 5 Page ! i f O ther Core M aterials ... . . . ..... . .. . . . . .. . . . . . . . . ... . .. . . .. . .. .... . . .. .. . .. ... .. ... .. . ... ...... 94 i L APPENDIX B-NCT GA....................................................................................95 ,
+
.i
)-
APPENDIX C - PARET... . .. .. . .... ....... . . . ... . .. . . .. . . .. .. ... . ... .. . . .. . .. . .. . ... . . ... . .. .. .. ...... .... . .. .... 99 , h A pPENDIX D . PROCESSING CODES ........ ................................................. ............ 112 i l l Power 2..................................................................................................I12 i AXPOW............................................................................................121 , GETFLUX............................................................................................126
- 7 VlTA....................................................................................................................129 ,
i 4 h l i i 1 i P I l l I; i s f I i
)
i t i s l I t l t L F h I I p i ? f I. I t I r i
1 l l , l vii l LIST OF FIGURES 1 l l ! FIGURE Page 1 l 1 N S CR gri d pl ate . .. . .. . .. . . . . . . . . . . . . . .. . .. . . . . . . . . .. . . . . . . . . . . ... . . . . . . . . .. . . . . 7 1 i ! 2 FLI P core configuration . ... . . . .. . . ..... ... . .. .. . .. . . .......... .. .... . . . ... . . ... .. 8 ! 3 Fueled followerinstallation.. . . . . ..................9 l l \ 4 Four rod fuel element assembly for the NSCR ..... ... ... ........ .. . . . . ... . . . . 10 I l 5 The instrumented fuel rod ............ ... . . . . . . . . . . . . . . . . . . . . . . . . . . .. . . . . . 12 l l i 6 Fueled follower fuel bundle assembly ...... .. ..... . .... ........ ..... .. . . ... 14 i 1 7 NSCR fueled follower control rod. .. ... . . . . . . . . . . . . . . . . . . . . . . . . . . . ... .. 16 1 8 FLIP pulsing characteristics measured from an experiment in 1993. . . .. ..17 ; 1 1 9 Prompt negative temperature coefficient for FLIP fuel .. ........ ... . ... ...... 19 1 10 Comparison of core powers from BOL and EOL $1.80 pulses 1 modeled in PARET .... .. ............. ... . . .....................................20 l 11 Comparison of fuel centerline temperatures from BOL and EOL
$1.80 pulses modeled in PARET...... .. . ...................................21 12 Nominal fuel rod spacing in the NSCR Core... . ....... ............ ... . . . . ... . . . . 2 3 13 Comparison of prompt negative tem, erature coefficients for GA test case. 28 14 NSCR core map profile in x-y plane as modeled in DIF3D .. .... .. .. . ........ 30 15 NSCR control rod positions in z-plane as modeled in DIF3D.... .. . ... ..... 31 16 Axial geometry for thermal hydraulic analysis......... . ....... ... .. ....... . . .. .. 33 1
1 I l 17 Axial regions for NCTRIGA ... . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ..... . . 3 4 l 18 Radial geometry for NCTRIGA. . . . . . . . . .. . .. .35 l I
viii Page 19 Radial fuel geometry used in PARET...... .. . ...... .......... ...............37 20 Axial nodes for PARET ............ . ..... . .... . .... ............................38 21 Integral temperature coefficients for FLIP fuel ........ .. . . . . . . . . . . . . . . . . . ... 40 22 Fuel centerline temperatures for the IF rod with varying flow rates calculat ed by PARET .... ... . .. . .... . . . .. . . . ... .. . . ... . .. .. .. . . . . . . .. ... . . . . . . . . ... .... . .. . . . . 4 2 23 Comparison of peak powers for pulses performed at the NSCR and modeled in PARET . ... ... ... . .. . .... . ... .. . .. ... . ..... .. . .. .. . . . .. . . . . ... . . .. . . . .. .. . .. . .. . . .. . . . . 4 3 24 Comparison of full width at half maximum of pulse peak for pulses performed at the NSCR and modeled in PARET .............. ...... ... .. ............ 43 25 Comparison of energy generated during pulses performed at the NS CR and mod eled in PARET.. ........ ............ .. ..................... ... .. ... .. ... . . 44 26 Comparison of peak pulsing temperatures from pulses performed at the NSCR and modeled in PARET . ..... . . ........... .... .............. . .. ... .. 45 27 Ratio of the peak fuel temperature calculated for the instrumented fuel rod using PARET to the peak thermocouple temperature measured for pulses performed at the NSCR . ....... . ..... . . . .... ............................. ..... 46 28 Power and energy from $1.65 pulse at 1400 MWD as calculated by P ARET . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .....................................47 29 Fuel centerline, fuel surface and clad temperatures for the instrumented channel from $1.65 pulse at 1400 MWD as calculated by PARET........... . 48 i 30 Relative power distribution of the NSCR at 1.49 MW calculated with D I F3 D . . . . . . .. . . . . . . . . . . . .. . . . .. . . . . . . . .. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 1 31 Relative power distribution of the NSCR at 1.0 MW calculated wi th DI F3 D . . . . . . . . . . . . .. . . . .. . . . . . . . . .. . . .. . . . . .. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 2
\
1 l Page 1 32 Axial distributions for the hottest elements at 1.0 MW and l 1.49 M W calcula ted with DI F3 D ........... ..................... .............. .............. .. 54 33 Fuel centerline, wall, and water temperatures for the " average" channel calculated with NCTRIGA at 1.0 MW and 1.49 MW.. . . .......... ... 55 l 1 34 Heat fluxes for the " average" channel calculated with NCTRIGA at i
- 1. 0 M W and 1. 4 9 M W . . . . . . . . . . . . ... . . . . . . . . . .. . . . . . . .. . . . . . . . . .. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 6 )
! 1 l 35 Heat transfer coefficients and Reynolds's numbers for the " average" channel calculated with NCTRIGA at 1.0 MW and 1.49 MW..... .............. 56 l l 36 NSCR core map with possible IF rod locations .............. ............... ............ 61 l l 1 37 Peak fuel temperature as a function of reactor power calculated i with N CTRI G A . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .. . . . . . . . . . . . . . . . . . 6 2 1 38 Relative power distribution calculated by DIF3D for the prepulsing case... 65 ) i
- 39 Power and energy calculated by PARET for $2.66 EOL pulse .............. ..... 66 40 Fuel centerline, fuel surface and clad temperatures calculated by PARET for $2.66 EOL pulse ................... ...... ....... ................................. ... 67 )
41 Peak fuel and thermocouple temperatures calculated with PARET for prompt reactivity insertions at EOL ....................... ................................ 68 l 42 Power and energy calculated by PARET for'an accidental pulse of
$ 1.96 from full power at EOL .............. .............................. .......................... 71 43 Fuel centerline, fuel surface and clad temperatures calculated by PARET for an accidental pulse of $1.96 from full power at EOL........... 72 44 Integral rod worth for the shim safety bank calculated with DIF3D..... . ... 74 45 Power and energy calculated by PARET for a control rod run out l from 1 Watt at EOL .... ...... ... .... ..... ... .. . .... .. . .. . . ... ..... .. ... ... ... . .. ....... ..... .. ..... 75 l
_ . _ . _ _ _ . _ . . . ~ _ . _ . _ . _ _ _ _ _ _ _ _ _ . . . . . _ _ _ . . . _ _ _ _ . _ _ _ _ . _ - i i X ! i Page , 46 Fuel centerline, fuel surface and clad temperatures for the hot channel { l calculated by PARET for a control rod run out from 1 Watt at EOL........... 76
)
( 47 Power calculated by PARET for a control rod run out j from 1.4 9 M W at EOL .. ..... . ... ... .. ........ ... . . .. ... . ........... ......... ... . ... .... .... .. . .. ... .. . 7 7 1 l 1 48 Fuel centerline, fuel surface and clad temperatures for the hot channnel l calculated by PARET for a control rod run out from 1.49 MW at EOL ...... 78 ! 49 Power calculated by PARET for control rod run out t- from 1.0 M W at EOL .... ... . . .... . . . . . . . . . . . .. ... . ... .. . . ... . . .... .. ... ..... .. ... .. ....... .. .. ... .. .. . 79 l y 50 Fuel centerline, fuel surface and clad temperatures for the hot channel i calculated by PARET for a control rod run out from 1.0 MW at EOL ........ 80 l i 51 TRIGA FLIP fuel rod geometry used in WIMSD4/m .................. . ............. 91 l 52 NSCR shim safety control rod geometry used in WIMSD4/m................ .... 92 l l l 53 Transient control rod geometry used in WIMSD4/m ................................... 93 i i 54 NSCR regulating control rod geometry used in WIMSD4/m........... ........... 94 ) 55 Comparison of axial distributions generated with AXPOW.......................125 i l l t I s o I
,m ..
xi LIST OF TABLES TABLE Page 1 Principal Fuel Element Design Parameters for TRIGA FLIP Fuel . ......... .11 2 Core Parameters from the GA Test Case........... ... ....... . ....... . ... . . ... .... 26 3 Peaking Factors from the GA Test Case.... . . .. . ... ................. . ...... ..... 27 4 Estimated Peak Thermal Flux at 2 MW 4-Rod Cluster TRIGA-LEU Reacto r . . . . . . . . . .. . .... . . . .. . . . . . . . .. .. . . . ..........................................28 5 Comparison of NCTRIGA Results to Experimental Data .. . . . ......... ......... 32 6 Integral Feedback Curve Fits ($) . ................................................40 7 Control Rod Worths Measured at the NSCR and Calculated with DIF3D... 57 8 Excess Reactivity at 300 Watts Measured at the NSCR and Calculated Using DIF3D... . . . .... . ... .. .............. .. 58 9 Shutdown Margin Measured at the NSCR and Calculated Using DIF3D.... 59 10 Temperatures and Powers Calculated Using PARET for Pulse I nsertions at EOL . . . . . . . . . . . . . . . . .. . . . .. . . . . . .. . .. . .. . . .. . . .. . . .. .. .... . . .. . . . . . .. . ..... .. . .. 67 11 Limiting Safety System Setting Values for Pulsing............ . ... . .... .... ... . 69 12 Summary of Radiation Exposure Following the Design Basis Accident .. .. 82 13 Atom Densities for the Fuel Rod Model....... ........... ............... ... .......... .... 90 l I 1
1 CIIAPTER I INTRODUCTION The Nuclear Science Center Reactor (NSCR) operated from 1962 until 1967 with MTR-type curved aluminum plate elements. During this time the reactor was operated extensively at a maximum power level of 100 kW. In 1968, the reactor was converted to General Atomic TRIGA type fuel rods and the maximum power level was increased to 1.0 MW. The initial core loading produced satisfactory operation, but the experimental capability was soon affected by both fuel burnup and samarium buildup. To restore excess reactivity, additional fuel was periodically added to the core and a graphite reflector was added to all four faces. This eventually led to a 126-element core with a resultant decrease in the flux of almost 40% and the elimination of most of the irradiation facilities. In August,1970, fuel followed control rods were installed to help solve the problem of maintaining excess reactivity. This installation required modification of the grid plete to allow passage of the fueled portion of the control rod through the grid plate. An average reactivity increase of $1.10 per fueled follower was achieved which extended the core life nearly two years. The high fuel burnup rate of standard TRIGA cores continued to be an operational problem for the NSCR which has operated at a burnup of approximately 100 MW-days per year since 1969. General Atomic FLIP (Fuel Life Improvement Program) fuel was designed to replace standard TRIGA fuel elements. FLIP fuel has the same dimensions as standard TRIGA fuel and differs only in material composition. The enrichment was increased from 20% to 70%, the hydrogen to zirconium ratio was decreased from approximately 1.7 to 1.6, and 1.5 weight percent natural erbium was added as a burnable poison. FLIP fuel has an expected lifetime of 9 MW-years which is considerably higher than the 1/2 MW-year life-This thesis follows the format of Nuclear Technology.
2 time of standard fuel. L l Since funds were not available for a complete FLIP core, the initial core operated with l a mixture of FLIP and standard TRIGA fuel. A precedent for this had been established by General Atomic when they operated a standard core loaded with eighteen centrally located Fl.IP elements in a fuel test program. Calculations were performed at Texas A&M which i led to the conclusion that satisfactory core arrangements were possible with a mixed core. This approach provided the additional advantage of achieving a substantially greater burnup in the standard fuel which was used during the operation of the mixed core. The I NSCR operated with two mixed core loadings containing 35 FLIP elements and 59 FLIP elements each between D73 and 1979. As funds became available, the amount of FLIP L fuel was increased until a complete FLIP loading was achieved. By 1979, all of the standard fuel had been replaced and Cores V-VIII were operated ! with FLIP fuel only.' From 1979 to the present, the maximum power level has remained at 1.0 MW. The capabilities of the NSCR can be enhanced by increasing the maximum ! reactor power to 1.49 MW. This is a feasible change since FLIP fuel has been licensed in t similar facilities for powers as high as 2.0 MW.2 To minimize the expense and complexity ; of this proposed power upgrade, the core configuration is expected to remained unchanged. 'I An increase in power level would create a higher neutron flux in the irradiation l facilities, and would allow for more versatile experimental operation. With a neutron flux of approximately 1.5 times the current level, sample exposure times could be reduced by one-third. The number of different types of experiments that the facility is capable of performing could also be increased, thereby generating more income to support reactor operation. As an example, a higher neutron flux in the neutron radiography facility would l result in higher resolution for the real time camera, and greater versatility in the types of objects that could be inspected. This thesis contains much of the information necessary for a safety analysis report i amendment to license a power upgrade. Information on core and fuel geometry such as the l
. . - ._. .. - . =_
f 1 3 1 arrangement of fuel rods, control elements and other components in the reactor, the composition of fuel and control elements are also included. In addition, the power densities, peaking factors, temperatures of the fuel elements, and the limiting conditions for operation are discussed. Physics parameters presented include the worths of each
- control rod, the shutdown margin, the maximum excess reac
- ivity, and the temperature coefficient of reactivity. Accident scenarios are also investigated. The physics parameters l are calculated numerically and verified for avuracy by comparison to standard test cases
! i or existing experimental data. The methods used to perform this analysis were obtained from work performed for the i study to change the core from the highly enriched FLIP fuel to low enriched uranium fuel.' l Several computer codes were used in performing this study. The first was WIMSD4/m, a one-dimensional neutron transport code capable of solving for fiux distributions in simple geometries using 69 fine energy groups. WIMSD4/m creates a spatially and energetically averaged set of homogenized cross sections in the seven standard energy groups used in TRIGA analyses for various subregions of the core. Since the fuel's neutron cross sections are highly temperature dependent, this code must be run for the entire operational temperature range during both normal and accident conditions. The seven group cross sections were used in the three dimensional neutron diffusion l code DIF3D. This code is capable of computing the kerr of the core for different l configurations, as well as the spatially dependent power and neutron fluxes for each energy group. Different control rod positions and fuel loading arrangements were investigated, and the prompt temperature coefficient of reactivity necessary for the transient calculations was generated as a function of temperature, j The licensing analysis for the power upgrade included a steady state thermal hydraulic study using the power distribution obtained from DIF3D. To perform the steady state analysis, the code NCTRIGA was employed. NCTRIGA is a one dimensional heat transfer code that calculates temperatures at specified axial positions along a single fuel
f B 4 rod charmel and determines the natural convection induced flow rate. This code uses the ; power distribution generated from the neutronic analysis, along with information on the geometry and material composition of the fuel, to produce temperatures in the fuel and to predict the coolant flow rate. The results from this analysis are crucial to the power upgrade, since the limiting safety facters are based on the temperature of the fuel at failure and departure from nucleate boiling at the surface of the fuel rod. j Studies of reactor transients at the upgraded power level needed to be performed to ; ensure safety during accident and pulsing conditions. The code selected for these analyses
~
was PARET, a multi-nodal point reactor kinetics n..o heat transfer code. PARET can be used to model reactor fuel and coolant response to control rod movement including the rod ejection excursions that need to be modeled to simulate reactor pulsing. This code was used to determine if off normal operations that were allowable at the old power level are also licensable at the upgraded power level. Using all of the tools mentioned above, along with a few processing codes that were written to facilitate repetitive processes and data transfer, many of the calculations necessary for the completion of a safety analysis report amendment were performed in this thesis. The validity of these codes for calculations of these types was verified by 1 comparisons with experimental data and/or accepted computational values, according to guidelines from the International Atomic Energy Agency and previous safety analysis reports.3 t l l
l 5 l CIIAPTER II ( REACTOR DESIGN l l The purpose of this thesis is to investigate the ability of the existing NSCR reactor to operate safely at a power level of 1.49 MW. Since this power upgrade will not involve any major core changes, a detailed description of the existing reactor systems is given in the rest of this chapter. In general, the NSCR operates with a core composed of 86 General Atomic FLIP fuel elements, four shim safety control rods with fueled followers, one regulating control rod, and one transient control rod. These elements attach to a reactor grid plate along with graphite reflecting elements, experimental notches and instrumentation. The entire core is covered by 26 feet of water. Detailed descriptions of these components are presented below. l Mechanical Design I The reactor core, control and instrumentation systems are supported by a bridge that spans the reactor pool. Mounted on four wheels, the bridge travels on rails at the sides of the pool. This allows the reactor to be moved from one operating position to another. Electrical power, control-circuit wiring, and compressed air are supplied to the core through the bridge. The reactor grid plate is welded to an aluminum suspension frame which is a welded structure of 3/8" x 2" x 2" aluminum angle. The west side of the frame is open toward the large section of the pool. The angle construction allows unrestricted flow of the cooling water. An aluminum stabilizer frame is bolted to the bottom of the grid plate for vertical support. Stainless steel guides on the bottom of the stabilizer fit between tracks on the pool floor. This allows accurate repositioning of the reactor core which is essential for l numerous experiments. The stabilizer also allows the core to be lowered to the bottom of i l the pool to prevent sway.
6 Fuel elements and control rods are contained in bundle assemblies which are positioned and supported by the grid plate in a 9 x 6 array of 54 holes. A top view of the grid plate is shown in Figure 1, The reactor core can be arranged in a number of ways on tne grid plate. A typical core loading containing 90 fuel elements and graphite reflectors is shown in Figure 2. In this configuration the A row of the grid plate is available for placing experiments. A 3/4" diameter clearance hole through the grid plate is used to allow passage of the fueled section of the control rods as shown in Figure 3. The NSCR currently utilizes General Atomic TRIGA FLIP type fuel moderator elements in which a zirconium hydride moderator is homogeneously combined with highly enriched uranium fuel. The expected fuel life of FLIP fuel has been stated to be 3000 4 MWD , although calculations support burnup as high as 4600 MWD. TRIGA FLIP fuel elements have a 15 inch long active fuel section. These elements contain approximately 235 8.5 weight-% uranium enriched to 70% in U that is homogeneous y mixed with ZrH and with approximately 1.5 weight-% erbium as a burnable poison. To facilitate hydriding, an 0.18 inch diameter hole is drilled through the center of the active section, and a zirconium rod is inserted in this hole after hydriding is complete. As shown in Figure 4,3.5 inch long graphite slugs, act as top and bottom reflectors. The active fuel section and top and bottom graphite slugs are contained in a 0.020-inch thick stainless steel clad. The stainless steel clad is welded to the top and bottom end fittings. 235 The approximate overall weight of the rod is 7 pounds with an average U content of about 123 grams. Serial numbers on the bottom end fittings are used to identify individual fuel rods. The principal design parameters are shown below in Table 1.
1 i l 7 N 1 1 A B C D E F o o l
, , l 9
i 1 l 8 7 I 6
- W 5 4
i 3 : i 2 , j . 1 O O , S ' Fig.1. NSCR grid plate. l
, 8 4 A B C D E F j 9
. o v o EM EGE i
. o : E R EN l OM :r g g :: ~. , p, sP QQ @ @t < f th O l
, , 2
^
, f- sp 1 .h @d POWER
& ^ @ t,1 3 g'A ,O 1 2, ,
E
.... g ,
j kg.: + Wk; J. q{,n w'yefg' v a 2 o n 1 - PN - - ; - a I i 1 l Shim Safety Rod With Fueled Follower h FlipFuel , h T ansient Rod With Air Follower h Instrumented Fuel h Regulating Rod With H2O Follower Pneumatic Tube E Experimenter Notch '*P ege tor Sb-Be Neutron Source Fig. 2. FLIP core configuration.
l ! 9
? '
rHOLD DOWN BAR
') UPPER GUIDE SECTION Ei b g s> '"
fFilELED FOLLOWER-
<i,' / POISONED SECTION *
$ln!
tk
/ ,
r
$iih
!!$ill.
E{:
> }
COOLING .: HOLES bg0 . i t
- LOWER GUIDE SECTION o
t{i k't i ' $!b CLEARANCE ( , j'*ill7 HOLE nbhha A i l f1 -#: ; i h, '17 i 4!!!!!!!!! i
~
R 1 .,. .ijjp *- ', en ; . n
*e ~
l ,
' r llll :
Mi e-- 4' niim-
"kfksiill!i jff.h [.:
i! ) gik, ih
$RF W FUELED FOLLOWER-FUELED SECTION l Fig. 3. Fueled follower installation.
l l 1 I
l l
~
'A % ,,,,,,,,,,, . . . -
l _;y
,. .)
_. . . . g. . _ _ __..__s g ... . ............ 4 , s YO O'
. ~
M-- cT1-1- a M
)b n
i M
; n; OYU - % ( L1_. '
.* ' J i
4 ROD FUEL ELEMENT ASSEMBLY n
) .kr b h TRIGA FUEL ROD ASSEMBLY t
Fig.4. Four rod fuel element assembly for the NSCR. C
i 11 TABLE 1 Principal Fuel Element Design Parameters for TRIGA FLIP Fuel Fuel-moderator material U-ZrII Uranium content 8.5 Wt-% U-235 enrichment 70 % U-235 content (avg.) per element 123 g i Bumable poison natural erbium Erbium content 1.5 wt-% ; Shape cylindrical i Length of fuel meat 15 in ) i Diameter of fuel meat 1.371 in i l Cladding material Type 304 SS Cladding thickness 0.020 in i l l 1 l A specially fabricated instrumented fuel element containing three thermocouples l embedded in the fuel is used to measure fuel temperature during reactor operation. As shown in Figure 5, tiic sensing tips of the instrumented fuel rod's thermocouples are i located halfway from the radial center line at the axial center of the fuel section and one inch above and below the axial center. The thermacuple leadout wires pass through a ] seal contained in a stainless steel tube welded to the upper end fixture. This tube projects t about three inches above the upper end of the element. A watertight conduit, composed of tubing connected by swagelock unions, guides the leadout wire above the water surface in the reactor pool. The NSCR fuel elements are arranged in two, three or four rod assemblies. The four rod assembly consists of an aluminum bottom adapter, four stainless steel clad TRIGA fuel rods, and an aluminum top fitting which serves as a handle. The bottom adapter fits into l
12 I t SOFT SOLDER STAINLESS STEEL
- FlLLER PLUG i
STAINLESS STEEL m i, j LEAD-0UT TUBE % ;!i' S U ;j i, SPACER END PLUG GRAPHlTE END REFLECTOR i ZlRCONIUM ROD 'j FUEL-MODERATOR I MATERIAL 0 r,. ( " T Ti 'h
- l THERMOCOUPLES (3) ,: : ..-N c
1 IN. 15 (I5,: l p l Q. IN. b:3. - ;
'1 IN.
4 O.02 IN. $.. N5, STAINLESS STEEL CLAOD1NG Mr.:k M 1 O. o GRAPHITE END REFLECTOR o l ! I i '. \ i 41 IN. Fig. 5. The instrumented fuel rod. r
l 13 l the NSCR grid plate, and the top adapter contains four tapped holes into which the fuel rods are threaded . The bottom fitting on the fuel rod has a flange at the base of the threads so that the fuel rod seats firmly on the adapter and is rigidly supported. The details of the four rod assembly are shown in Figure 4. A three rod fuel assembly may accommodate a control rod, an instrumented fuel element, or an experiment. The NSCR utilizes two separate types of three rod fuel . assemblies for housing control rods. The first permits one fuel rod in an assembly to be t , ! replaced by a control rod guide tube which has an outside diameter of 1.5 inches. The ! handles on control rod elements are modified to accommodate the guide tube. A l regulating rod and a transient rod without a follower will utilize this type of assembly. The second type, shown in Figure 6, is designed for use with shim safety control rods which are fuel followed. The transient rod which has a follower uses a specially designed control rod guide tube and must also have a base assembiv as shown in Figure 6. The instrumented fuel rod fits into the bottom adapter. Not being an integral part of the bundle, the instrumented fuel rod is positioned into the bundle after it is in the grid plate. Two rod assemblies are identical to four rod assemblies except that they contain only ' two fuel rods. These assemblies are used to produce' thermal neutron flux traps by replacing fuel with water. Two rod bundles are positioned along the core's periphery in order to enhance the neutron flux available at certain irradiation locations. The nonfuel section of a two rod assembly may also be used to contain an experiment. Six motor-driven control rods (four shim-safety rods, a regulating rod, and a transient rod control the reactor during steady-state operation. Typical control rod locations are shown in Figure 2. The shim-safety control rods have scram capability and contain borated graphite powder as a poison in a stainless steel cladding. The shim safety control i ! rods consist of a 35.56 cm (14 in.) long poison section with a diameter of 1.7 cm en top of a 38.10 cm (15 in.) fueled section with a composition identical to the fuel rods, but with a fuel diameter of only 1.7 cm instead of the 1.74 cm fuel rod diameter. These are capped l ! 1 l
14 l FUEL ELEMENT HOLD-DOWN BLADE LOCKING BOLT HOLD-DOWN FOOT
/k
/f N o C o o
LOCKING PLATE l ~o UPPER GUIDE i n
.Q Q i
d d FUEL ROD , T i l @v C gFUELED FOLLOWER O /
/
l~f l~{ E E v - x - LOW ER GUIDE T
% s
% mi Fig. 6. Fueled follower fuel bundle assembly.
i
15 on both ends by void regions as shown in Figure 7. The regulating rod consists of a 38.10 cm (15 in.) poison section with a diameter of 1.52 cm. The transient rod has a 38.10 cm (15 in.) poison region on top of a 50.80 cm (20 in.) void section. Both of these sections have diameters of 1.52 cm. The shim safety rods and the transient rod have scram capability, but the regulating rod does not. The shim safety and regulating control rods are moved axially throughout the core with electromechanical drive mechanisms. The shim safety rod drives are controlled manually, and can be moved individually or in gang. Their maximum withdrawal speed is 11.4 cm/ min. The regulating rod can be controlled manually, or operated in automatic mode, during which the rod position is varied to maintain a constant power level. Its maximum withdrawal rate is 24.4 cm/ min. The transient rod is held in position by high pressure air, and is capable of rapid ejection from the core for pulsing. Nuclear Design The design and operating characteristics of standard TRIGA cores are well known as is the inherent safety characteristics of this class of reactors.5 TRIGA fuel is designed with physical parameters to insure safety under a wide range of operating conditions. Pulsing is one of this reactor's most notable capabilities. When the transient rod is ejected, the reactor can reach powers as high as 1000 MW for millisecond time intervals; as heat begins to build up in the fuel, the power automatically decreases due to the prompt negative temperature coefficient of TRIGA fuel. The current reactivity limit for pulsing the NSCR is $2.35 which results in a safe, licensable condition. It is worth noting that earlier operational NSCR cores were regularly pulsed with $3.00 insertions. The pulsing characteristics for an operational NSCR TRIGA-FLIP core are shown in Figure 8 for pulses as high as $1.80. This data was l , recorded during an experiment performed in 1993. l i i
1 i i If 1 l 1 4 l li "ft) I , l 9
/-V'D ;
1 ! l l 634 7 375S.S. 0.D. x .02 WALL TUBING I rg S.S. PLUG f fBORATED GRAPHITE 14 { - y ' 1 g i, d!I l 3g , S.S. PLUG
,j 7 ZlRCONIUM ROD 4 <
g FUEL l 45j
~
- 15 l
l r S.S. PLUG l g VOID , 1
]
7j : hN l Fig. 7. NSCR fueled follower control rod. i
_ . _ . _ _ _ . . . . . . _ _ . . _ . _ . _ _ . _ . . _ ~ _ _ = . . _ . _ _ _ _ - . . _ _ . . _ _ . _ _ _ . _ . _ __ _ 17 ; i l 1000 m i
<T !
../ i G /' ^ b
? - ^ /
~
?
3 h.-- / , h C OPeak Power P 100 , 3-GEnergy ; g t-r 7' A Apeak 1hermocouple Temp. A N m- M] bU JW
~
X 0 10 m' > F : E t
.g 1
1.0 1.2 1.4 1.6 1.8 j Reactivity Insertion ($) Fig. 8. FLIP pulsing characteristics measured from an experiment in 1993. TRIGA FLIP fuel contains two important constituents that allow pulsing to occur ; safely. One of these is erbium with natural isotopics which is loaded to 1.5 weight % as a burnable absorber. i67Er has a large resonance cross section at ~0.5 eV. As the neutron spectrum is hardened by fuel and moderator temperature increases during the pulse, more 167 l neutrons reach the energy range of the Er resonance and are absorbed. TRIGA FLIP fuel also contains ZrH as a moderator in the fuel structure which allows fission neutrons to be thermalized in the fuel rod where they are produced. As the fuel temperature increases, the ZrH heats immediately and moderates the neutrons to a higher average energy, placing i 67 Er resonance and resulting in fewer fissions. The ZrH and even more neutrons in the
'67
- Er together enhance the negative temperature feedback.
f
is The calculated temperature coefficients for FLIP fuel are shown in Figure 9. The temperature dependent variation of the temperature coefficient of a TRIGA FLIP core is ! advantageous in that an acceptable reactivity loss is incurred in reaching normal operating l temperatures, but any further increases in the average core temperature result in a sizably increased prompt negative temperature coefficient to act as a shutdown mechanism.
-l 235 Burnup calculations indicate that after 3000 MWD of operation with FLIP fuel, the U !
i 167 concentration averaged over the core is -67% and the Er concentration is ~33% of the { i beginning-of-life values. The temperature coefficient of reactivity at end-of-life is less ! temperature dependent than the beginning-of-life temperature coefficient of reactivity j 367 because of the sizable loss of Er and the resulting loss of the ~0.5 eV resonance absorption. l ! The transient behavior during pulses of $1.80 are shown in Figure 10 for both ! beginning oflife and end oflife conditions. The effects of the decrease in the erbium ' concentration can be seen since the same insertion results in a higher power when the pulse ! l occurs at the end of life. The peak fuel temperatures for these pulses are shown in i
-i Figure 11. As expected, the peak temperature at end of life is higher than that'at the j beginning-of-life. !
i Fuel Description and Safety Limits ! The most restrictive limit on the operation of a TRIGA core results from the l dissociation of hydrogen from the uranium-zirconium hydride fuel. As the hydrogen is f dissociated from the fuelmatrix, pressure increases inside the fuel rod and can eventually j i- ! cause a cladding failure. Another concern is the migration of hydrogen from the hot radial ! center of the fuel to the pellet edge under prolonged irradiation. Hydrogen to zirconium - ! ratios near the edge of the fuel can increase by over 15% from their initial values. During , i .
- a pulse, the temperature is momentarily peaked at the edge of the fuel, and even greater L
1 l i i
1 1 I l 19 l 1 4 I i l l l 18 0 MWD i
- - - -1000 MWD 1 16
- - - - - - 1400 MWD ;
_ . _ . -2000 MWD 1
- - - - -3000 fAND j 1
14 i
~
b. O 12 - 1 h . l.i .
'c 10 E -
g , ~ ~ _ _ , ; U 2 ti; 8 / e ss' *-- y
/
/
/
,,'"..**,,-- "' ~.,,
u .* k6 [, * #
/ .', s. '. _ ~
--~-
~ ~~
1 I i "" .'* ' E
/ - %
l-4
/,'dl//,,// " .%*N N. '
/ ~%.. %
~.
2 % I - 0 0 100 200 300 400 500 000 700 800 Fuel Temperature ('C) Fig. 9. Prompt negative temperature coefficient for FLIP fuel. ; 1 i e I i 1 i
- l. 1
i i 1000 . . . . > > i 900 - ,5 i'
'tl 800 -
l1' - 700 -
'l
'h I
600 - k2 i
' BOL -
500 -
- - - - EOL f,
O ' A - 400 - i i I 300 - l 1 I 200 - l i 100 - i i 1
- - - ~ ~ - - - ~ ~ - - - - ~ - - - - " - - - - ' - - ^ - - - ~
0 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 Time (seconds) Fig.10. Comparison of core powers from BOL and EOL $1.80 pulses modeled in PARET. E$ oee n w er. 3+-w- .--_---.--.m --y,w.. e- -- e r.*>e-- .e-e. e- - -- p.- -- -- eg m o. .2 h
700 . . . . . , . . 650 -
~~~~~~~~~~~~ -
600 - -
/
_- . . ~~ ' 550 - ,
~~ %
/ u 500 - ,
~
, 's -
G 450 -l , w . _ y 400 ,' o - s 350 -l y E ' 300 -j BOL
- - - EOL _
# 250 _i 1,
200 -t I - 150 i 100 -! , 50 -l 0 12 13 14 15 O 1 2 3 4 5 6 7 8 9 10 11 Time (seconds) Fig.11. Comparison of fuel centerline temperatures from BOL and EOL $1.80 pulses S modeled in PARET.
22 hydrogen dissociation can occur. The resulting pressure can cause b element deformation. Pulse insertion limits may have to be adjusted at high prevent this deformation. Studies show that hydrogen pressures fits ultimate temperature of up to 1150 C would not produce a stress in the clad in excess In order to reduce strength as long as the cladding temperature does not exceed 500 C.6 limit in the risk of fuel deformation, all operations must be such that the peak temperatu any fuel rod, under all conditions of operation, will not to exceed 950 C. Thermal Design The NSCR is cooled by natural convection and can be operated at a i The the pool center line except near the gateway between theh stall and lar reactor core is constantly surrounded by pool water which is drawn in f ; bottom and sides of the core during the convective cooling process. ' Figure 4 shows that the four rod assembly has been designed to pro l Water is drawn into the assembly by natural of cooling water through the element. It passes through the convection through the 2" diameter hole in the grid plate adapter. large cruciform opening and then over the entire element until the numerous openings in the aluminum handle. In addition to the coolan I h diameter l through the grid plate adapters, the NSCR grid plate has additional 1 coolant holes located at the corner of each four rod element. l Mark 111 standard fuel elements and FLIP clements have been TRIGA cores by General Atomic at steady state power ig levels of up i l arrangement of fuel in the NSCR has been designed so that the minim between the fuel rods provides adequate convective cooling for The spacing of the fuel rods in the NSCR core is shown in F be enhanced by the rod to rod spacing and by the extra cooling h
! l l
23 l I bundle. Cooling of the NSCR is also improved due to the increased depth of pool. The j NSCR core is normally covered by 26 feet of water which places it at a depth greater than 4 i that of most TRIGA installations. The increase in hydrostatic pressure created by this , I j greater depth reduces nucleate bubble formation and improves the margin to departure to - j nucleate boiling (DNB). 1
& 3.189 h. 8 4
i i e4g g> +e. <
. . 1330in.
l
. . 0094 in.
[ . \! - V 4 . .
-{- . 1.530 n
~ ~
\
.530 in.
GRID ILATE P G ADAN 0248 in. ! l Fig.12. Nominal fuel rod spacing in the NSCR core. I i l l l
i 24 1 CIIAPTER III NUMERICAL ANALYSIS TECIINIQUES The investigation of a core composed of FLIP fuel required the use of computer codes capable of producing all of the necessary operational parameters during both steady state operation and accident transients. The accuracy of computations is of the utmost importance, and is assured by performing benchmark calculations. This chapter contains a description of the methods used to produce the results presented later in this thesis along with a verification of these methods. Neutronic Analysis The neutronic analysis performed in this thesis was based on techniques developed and described in the report on the conversion of the NSCR from FLIP to LEU type TRIGA fuel. ! Many of the same cross section files and core configurations used for the LEU conversion analyses were also used in this thesis. Although the techniques were not developed for exclusive use here, they are presented in detail below for completeness. The code selected for the generaticn of few group neutron cross sections representative of each subregion in the core and its surroundings is WIMSD4/m, a one-dimensional neutron transport code capable of solving for flux distributions using 69 fine energy groups7 . WIMSD4/m is used to create spatially and energetically averaged sets of cross sections in the seven standard energy groups used in TRIGA analyses for each subregion of the core as a function of fuel burnup. Since the neutron cross sections in the fuel are highly temperature dependent, WIMSD4/m was run for many cases covering the entire operational temperature range for both normal and accident conditions for each type of subregion in tne core. WIMSD4/m is especially suited for the analysis of TRIGA fuel since it can model neutron scattering with hydrogen in ZrH. It's isotope library also 166 contains cross sections for Er and 67Er.8 Since extensive sets of cross sections for TRIGA FLIP fuel and the other core materials for the NSCR had been previously l
\
25 generated in the LEU conversion report, the details of their generation are not presented here, but are included in Appendix A. Once a temperature and burnup dependent library of homogenized cross sections was generated, it was used in a three dimensional neutron diffusion calculational model of the core, reflector and irradiation facilities. The code selected for these calculations was DIF3D'. This code is capable of computing the kg of the core for different fuel and control rod configurations, as well as the spatially dependent power and neutron flux distributions for each energy group. Different control rod positions and fuel loading arrangements can be investigated, and the prompt temperature coefficients of reactivity necessary for the transient calculations can be generated on a full core calculational basis. To prove the accuracy of these codes for modeling the NSCR, a test case presented by General Atomic for a two megawatt TRIGA core was recreated.3 Although the fuel modeled in this test case is specifically TRIGA-LEU fuel which contains 20 weight-% 235 uranium enriched to 20 weight-% in U, many of the properties of the fuel are similar to TRIGA FLIP fuel, including the physical dimensions, the use of Zrl! as a moderator in the fuel, and the use of Er as a burnable poison. This test case has been accepted by the International Atomic Energy Agency to be representative of the properties of TRIGA cores.3 This problem was chosen because it has clearly defined conditions and detailed l results to allow verification of many important neutronic parameters. A summary of the core parameters for this test case is presented in Table 2.
p 26 I TABLE 2 Core Parameters from the GA Test Case l i Fuel cluster: TRIGA-LEU 20 wt-% U in UZrH (76 x 80 x 508 mm)
- Fuel rods per cluster
Standard Cluster: '4 ( Control Cluster: 3 Nominal fuel rod dimensions: Fuel O.D.: 32.4 mm Clad O.D.: 33.5 mm (incoloy) Fuel Height: 508 mm
- Fuelloading: 548 mm U (20% enriched)/ rod 2.2 Kg U (20% enriched)/std cluster 440 gm U-235/std cluster 4.59 wt-% Erbium as burnable absorber Number of fuel clusters in the core: 26il Standard clusters: 21 Control clusters: 5il Refhetor: Water Core size (liters): 78i2
~ U-235 content / core (Kg): 10.6 Core geometry: 4 x 6 arrangement Grid plate: 6 x 9 positions (normal conversion)
Bumup Status of the core: equilibrium core Average core burnup (%): ~20 . Thermal-hydraulic data: Average power density (kW/ liter): 26 - Coolant flow rate: 1000 gpm (227 m3 /hr) Core inlet temperature: 38 C l
27 A set of peaking factors were generated and reported by GA for the test case. These are listed below in Table 3 along with the corresponding values calculated in this thesis , using WIMSD4/m and DIF3D. TABLE 3 Peaking Factors from the GA Test Case Type of Peaking P/P General P/P TAMU P/P TAMU Atomic 0 MWD 900 MWD Core Radial 1.57 1.61 1.62 Core Axial 1.36 1.30 1.31 1D Cell (23 C) 1.48 1.41 1.29 1D Cell (310 C) 1.52 1.44 1.31 1D Cell (700 C) 1.61 1.49 1.34 The core radial and core axial peaking factor values agree quite well, but there is a l small discrepancy in the ID cell values. The GA test case was run at a burnup of approximately 20% of the total core lifetime which corresponds to approximately 900 MWD. It is possible, however, that the ID cell peaking factors presented in the document were calculated with no bumup on the core. This would put them in better agreement with i I the TAMU values. Table 4 shows a comparison of the peak thermal flux values in the core and the water reflector. To find the thermal flux, the values of the flux in groups 5,6, and 7 in the seven
- group model were summed. This includes all neutrons with energies below 0.5 eV. The l
resulting flux values can be seen to be in close agreement. t
28 TABLE 4 Estimated Peak Thermal Flux at 2 MW 4. Rod Cluster TRIGA. LEU Reactor General Atomic TAMU 3 Core 1.5 x 10'3 1.5 x 10 3 Reflector (water) 2 x 10 3 1.9 x 10 The final parameter that GA reported for the test case was the prompt negative temperature coefficient. Figure 13 shows the TAMU values along with the GA values for this parameter. p 18 2 - y 16 - ' h .-
$ 14 ,'
E
$ 12 ,.
Q
- j 10 ,.
b , TAMU
'f 8
,',' ......og N ,
y6-8 g4 b
$2-e 0 l l l -
0 200 400 600 800 1000 Temperature (deg C) l Fig.13. Comparison of prompt negative temperature coefficients for GA test case. l l 1 j
i , 29 l Although these are the only parameters which could be verified against General Atomic results for the benchmark case, they show that the TAMU models are capable of
- producing accurate results as compared to the GA approved models. DIF3D is also l
capable of producing output in much more detail than that given here. With the aid of an external plotting program (TECPLOT), three dimensional flux maps can be generated which are useful for visualizing the neutron population in the core. l The DIF3D model of the NSCR consists of parallelepiped blocks of the various j materials with cross sections that are defined in Appendix A. A diagram of the material arrangement in one K-y plane at the axial centerline is shown in Figure 14. This arrangement is surrounded on all sides by water. The contents of the planes in the z-direction of the three dimensional model can be varied to model changing material configurations, such as different control rod positions. Some sample axial configurations are shown in Figure 15 for FLIP fuel and control rod positions. Control rod movement can only be modeled in steps of one inch each; in restity, the control rods move continuously. To verify the above DIF3D model of the NSCR, a test case was run and compared to measured critical data taken in December 1993. The core average temperature and bumup were assumed to be 500 K and 1400 MWDs respectively, and the control rod positions were as follows: shim safety rods at 66% out, transient rod at 100% out, and regulating rod at 45% out. DIF3D converged to a kg of 1.011 for this case. This gives a bias of $1.57, 235 but since the choice among the various U cross section sets in the WIMSD4/m library change the calculated ken, the fuel rod burnups could not be verified, and the control rod ! modeling was approximate, a ker 1.011 has been accepted as " critical" for this model. All 1 other cases will be considered critical when DIF3D converges to a kwro f 1.011. ! l l
J i 1 j 30 i i i i i 1 .i l E I i 1 Shim Safety Rod No.1 i 16.192cm i 1
][
,i _2_ Shim Safety Rod No. 2 j h i ,
. i.
.. 3 Sli. dafety Rod No. 3 l ;: ; "b 3, . \ . : ;}. -
l ! 3; i .. 4, , .1 _4_ Shim Safety Rod No. 4 i 7 l-40.480cm - - 1 bbb...; . L IL..;L. i w L-.. TR Transient Rod
..7 L_.A_ ;L2 L_. L.i, . . t .. . u c. _:
.I.'..; RR Regulating Rod
. . - .:.; 1 . _ A; r*- ,
p V bbbb [ .b kI..} j- 7 Fuel l , . Graphite l r f 4 16.192cm 1 - 1 p Water s, f- s_ . MM 38.560cm M l 7.712cm 4 A ! Fig.14. NSCR core map profile in x-y plane as modeled in DIF3D. i i $, 1 i l 1 1 J i i 4 4
i 31 a 1 i i i I Rods Fully Inserted Rods Fully Withdrawn Up JL L i 22.86cm JL i ll.43cm 2b .- ,,,,,, , ,
,,,, ,,,, ,,m 4 N 8.89cm ; l U
i JL I E ControlRodPoison ,
] ControlRod Void l l ,
g Fuel y . : 8.89cm 5
- g Graphite End Reflector ,
p of Fuel Rod j i, , , , ,, , , ,i,, 8.89cm @ Top and Bottom j j[ l Stmetural Material JL f 12.70cm g - g Aluminum Grid Plate ! ! JL 12.70cm i O w ater 1 3 efy Trans. Reg. 3 ey '* Trans. Rod Rod Rod Rod Rod RRed i Fig.15. NSCR control rod positions in z-plane as modeled in DIF3D. I
- l i
e
l 32 ThermalIlydraulie Analysis A steady state thermal hydraulic analysis must be performed using the power distributions obtained from DIF3D to determine the maximum fuel and clad temperatures during operation. To perform this steady state analysis, the code NCTRIGA was used. NCTRIGA is a one dimensional thermal hydraulic code that calculates the temperatures of the fuel centerline, fuel surface, and coolant at several nodes along a single fuel rod channel. It also determines the natural convection induced coolant flow rate. This code uses the power shape generated from the neutronic analysis, along with information on the geometry and material composition of the calculational channel, to calculate temperatures in the fuel and to predict the coolant flow rate. This code is based on the Argonne National Laboratory's NATCON code and has been modified to model reactors with TRIGA type fuel. ' Verification of this code was performed at Argonne National Laboratory with data collected by General Atomic. The data shown in Table 5 verifies that reasonable temperatures and flow rates can be expected when applying NCTRIGA to TRIGA cores. TABLE 5 Comparison of NCTRIGA Results to Experimental Data Power No. of Flow rate Error Outlet Error MW Elements Source kg/sec % Temp C % 1.0 91 GA 8457 70.2 NCTRIGA 7890 -7.0 67.2 -10.0 l 1.5 91 GA 9555 76.6 ; NCTRIGA 9259 -3.2 75.6 -2.6 j 2.0 101 GA 11080 86.1 ! NCTRIGA 11051 -0.3 80.1 -14.0
- NCTRIGA will be used to find all of the steady state temperatures needed to produce l
j the limiting safety system settings and to find the initial flow rates necessary for the transient analysis. The input for NCTRIGA is relatively short since only one fuel rod is j I
33 modeled. Detailed information on the whole core is not necessary; only the total core power and number of fuel elements are needed. NCTRIGA requires basic information about the fuel geometry and material properties, along with the axial power distribution and radial peaking factor for the rod ofinterest. The geometry of a " channel" is modeled in NCTRIGA as a fueled region capped on the top and bottom by non-fueled regions. The NCTRIGA geometry is shown below in Figure 16. A detailed description of the input parameters and a sample input deck are given in Appendix B. x-3 2032 m i 45 l Ftel Region 5 Non-fuelRegion I 3 381 m iH 31778 m
-m Fig.16. Axial geometry for thermal hydraulic analysis.
I i
l l 34 Since the diffusion code generates power distribution data for each of 15 axial regions, the same regions were used in NCTRIGA. NCTRIGA requires the power distribution to be input for the ends of each axial region as shown in Figure 17 below. 16 Is 15 j 14 I4 13 Geometrical 13 Power 37 Regions 12 Distribution 11[ 10 Nodes 9 R 8 7 6 5 5 4 4 1 2 2 1 1 Fig.17. Axial regions for NCTRIGA. To generate the axial power values at these nodes, a processing code AXPOW was written. AXPOW linearly extrapolates the DIF3D data to generate the power values at the top and bottom of the fuel region, and perfomis a higher order interpolation using Neville's Algorithm for the central nodes. AXPOW then prints the data in the format necessary for input to NCTRIGA. A source listing for AXPOW is given in Appendix D. Radially, a " channel" is modeled in NCTRIGA as a cylindrical fuel rod with a radius of 0.017411 m, surrounded by a cladding with a thickness of 0.000508 m. The fuel clad gap is represented in NCTRIGA by its conductance and is not given a physical thickness. l
l 35 . i l The pitch to diameter ratio of 1.10268 gives the necessary information for the calculation of the water channel dimensions. The radial geometry is shown in Figure 18. l Fuel l j _-~_
~
Clad ! 1 E l 1
^
l Fig. I8. Radial geometry for NCTRIGA. To perform calculations to determine the peak fuel and clad temperatures during pulsing and accident transients, the ANL's PARET code was used. PARET is designed for use in predicting the course and consequences of nondestructive reactivity accidents in small reactor cores. It is a coupled neutronic-hydrodynamic-heat transfer code employing point kinetics, one-dimensional hydrodynamics, and one-dimensional heat transfer." The core can be represented by one to four channels. Each of the four channels may have different power generation, coolant mass flow rate, and hydraulic parameters as represented by a single fuel pin and its associated coolant flow. The time dependent temperatures within the fuel element are computed on the basis of the one-dimensional heat conduction equation solved in each of up to 21 axial sections. The hydrodynamics solution is also one-dimensional for each channel at each time node. The heat transfer may ! take place by natural or forced convection with nucleate, transition, or stable film boiling, and the coolant is allowed to range from subcooled liquid, through the two-phase regime, up to and including super-heated steam. PARET will also calculate flow reversal.
36 To model the NSCR, a PARET deck was set up with two representative channels. The first is the channel ofinterest, referred to as the hot channel. It represents the channel with the highest radial peaking factor or the channel containing the instrumented fuel element. The axial power distribution for this channel was taken from DIF3D for the core location ofinterest. In PARET, the hot channel is given a reactivity feedback weighting factor of 1/90, so its temperature will contribute 1/90 of the core average temperature in the reactivity feedback calculations. The second channei represents an average fuel rod in the core. The axial power distribution for this channel is taken from the DIF3D results for a channel with a radial peaking factor close to one. The reactivity feedback weighting factor assigned to this channel in PARET is 89/90, so that this channel represents the remaining fuel in the core during the temperature feedback calculations. A detailed description of PARET input can be found in Appendix C. A summary of the primary parameters is given below. The geometry of each channel is modeled in PARET as follows. Axially, the channel is divided into three regions. There is one fuel region, and two non-fuel entrance and exit regions which represent the graphite reflectors, and end plugs on each fuel rod. The axial geometry is shown in Figure 16. The entrance and exit regions of the channel are represented in PARET by their equivalent hydraulic diameter. The outer radius of the fuel rod is 0.01792 m over the axial length of the channel, including the non-fuel regions. The fuel region is divided into six radial regions. This division followed the recommendations in the PARET manual. The first four represent the U-ZrH fuel. The fifth region represents the gap, and the last represents the cladding. The fuel regions are each 0.00435275 m thick. The gap is 9.0x10 4 m thick, and the clad is 5.0x10" m thick. These are surrounded by a coolant region which has a thickness of - l 0.01 m. This is shown below in Figure 19. WIMSD4/m was used to generate the radial within rod power shape for PARET. This - power shape was generated from the thermal flux edit of WIMSD4/m. The code
l l ; l ! 37 l GETFLUX reads the WIMSD4/m output file, interpolates the data, and performs volume weighted averaging of the thermal flux over the cells used in PARET. This code is showTi l l in Appendix D. i l In order to input the axial power shape for each channel, the total fueled region must be i subdivided into smaller axial regions. In this model, each region is taken to be one inch in length, or 0.0254 m, to correspond to the 15 axial regions used in the DIF3D model. PARET requires that the power shape be given at the ends of the top and bottom regions l : I and at the center of all of the other regions. This is shown graphically in Figure 20. Ash 2 l Th Fuel e ',
)
7 ,' ,
, yGap
[ Alad
'h j\ , >
[' i p q - - - j oolant Fig.19. Radial fuel geometry used in PARET. 1 l l l i
i l 38 15 15 14 14 13 13 Geometrical Power D 12 Regions Distribution ! f 33 Nodes ! 10 10 9 9 l 8 8 7 7 6 6 5 5 : l 4 4 3 3 2 2 1 1 Fig. 20. Axial nodes for PARET. t Again, the code AXPOW is used to produce the power values at each location. The i values at the top and bottom nodes are produced by linearly extrapolating from the DIF3D l data, and the values at the central nodes are read directly from the DIF3D output. Since PARET has no separate input for peaking factors, all values must be multiplied by the l radial peaking factor for the channel ofinterest. AXPOW performs these calculations, and creates output in the format necessary for input to PARET. The fuel properties needed by PARET are the thermal conductivity, heat capacity, and I temperature reactivity feedback. The thermal conductivity of the fuel meat was taken to be constant with respect to temperature. Its value is given by GA as 18.0 bnK.i2 The thermal conductivities of the gap and clad are taken as 0.199 W/mK and 16.8 W/mK, respectively. The heat capacity is calculated as follows using techniques suggested by GA. The specific heat ofZrH i.6 is given by Cp(ZrH) = 0.75148 T + 363.0938 J/kg*C (1) l 4
l ! 39 1 The specific heat of U is given by Cp(U) = 0.1305 T + 109.4 J/kg C (2) The specific heat of UZrH is then found from Cp (UZrH) = w(U)Cp (U)+ w(ZrH)Cp (ZrH) (3) l where w(x) represents the weight percent of constituent x. Combining equations (1), (2), and (3) l Cp(UZrH) = w(U)(0.1305 T + 109.4) + w(ZrH)(0.75148 T + 363.0938) J/kg C (4) l Applying the appropriate weight percents of U and ZrH for FLIP fuel and multiplying by the density, 3 Cp(UZrH)=((0.085)(0.1305+109.4)+(0.915)(0.75148T+363.0938)J/kg C)x5970 kg/m (5)
~
finally gives il1e following result: 4 6 3 Cp (UZrH) = 0.417x10 T + 2.05x10 J/m C for FLIP fuel (6) The temperature coefficients of reactivity needed by PARET were generated using DIF3D with temperature dependent cross sections taken from WIMSD4/m. These are shown in Figure 21 for FLIP fuel at various values of burnup. 1 I l
i h i ! t i 40 t I 14 - - - - - - - - - - - - 0 MWD ,
- - - - 1000 MWD ;
6 12 - - - - - - 1400 MWD ! g - - - - -2000 MWD , l 2 - - - - -3000 MWD i i j 10 - ! w I 8- i E. - :
$ /
j e- f. ' i
? /y/
y*' //, . ', s # h 2 0 _fh I O 200 400 000 800 1000 1200
- 1 Temperature (K) :
t Fig. 21. Integral temperature coefficients for FLIP fuel. : To prepare the data for PARET input, curve fits were first generated. It was found that I second order fits gave an acceptably small error. These fits are shown below in Table 6 for i various core average burnups. ! s ! TABLE 6 Integral Feedback Curve Fits ($) l 4 2 BOL F(T) = -0.99258 - 7.57122x10 T + 1.26832x10 5 T 1000 MWD 4 2 F(1) = -1.63758 + 3.40290x10 3 T + 5.25240x10 T 2 1400 MWD F(T) = -1.92043 + 4.96090x10-3 T + 2.90688x10-6 T 2 2000 MWD F(T) = -2.27714 + 6.87820x10-3 T + 1.30048x10'7 T 3000 MWD 4 2 F(T) = -2.66930 + 8.86260x10-3 T - 2.43340x10 T
41 To acquire power distribution data for the pulsing model, a DIF3D job was run at 300 watts steady state power with the transient rod fully inserted to obtain the initial power distribution. This data was compiled into a PARET input deck to produce a model for predicting the performance of the peak fuel element and " average" fuel element during pulsing. The mass flow rates needed for input to PARET were generated using NCTRIGA. The same initial power distribution was input to NCTRIGA which was executed to produce the mass flow rate. This value was then put into the PARET input deck as the initial flow rate. Su'osequent values of the mass flow rates have little effect on the temperatures generated during the short pulsing transients, as long as sufficient flow is provided to allow the code to run to completion. Several jobs were run to verify this. The fuel centerline temperatures for the IF rod channel from these cases are shown in Figure 22. To benchmark the PAIET model, data was obtained from pulsing experiments on the NSCR core using FLIP fuel. These experiments were performed on April 8,1993. The ; average core burnup was approximately 1400 MWD. DIF3D was run to create a critical case with the same core configuration as the experiment. The power distribution of the i f instrumented element and average element were input to the PAIET model. The curve fit for the temperature feedback at 1400 MWD, as shown in Table 6, was also used. r l l I I l
i 600.0 i . i 1 e i G 400.0 - - U { 500 kg/m2s _ti 1000 kg/m2s 2 - - 1500 kg/m2s 8. e a e il a 200.0 - - i i 3 0.0 O.0 5.0 - 10.0 15.0-Time (seconds) Fig. 22. Fuel centerline temperatures for the IF rod with varying flow rates as calculated by PARET. _8 __._______._______________________________._aw-, m _m_m_-_ _..+wwea e -e. _ 10 _ g-PARET i .05 0 + 1.03 1.10 120 1.30 1A0 1.50 1.60 170 120 i Reactivity Inserted ($) Fig. 25. Comparison of energy generated during pulses performed at the NSCR and i modeled in PARET. The low value of the experimental data for the $1.80 pulse was most likely caused by ; l the detector saturating due to the excessive amount of energy generated in this pulse, and it ; r is therefore believed to be incorrect. The final comparison of the PARET model and experimental results considered the l peak temperature generated in the fuel during a pulse. These results are shown below in Figure 26. ! 4 t 1 1 1 ! i
45 1200 m .g U ' 400 O l M TEei~ M F ~" 0 -E-PARET cent line
-4 _PARET fuel surfa j 0- - -+- 4 * - - - -- +----4---
1 00 1.10 1.20 1 30 1 40 1.50 1.60 1 70 1 30 Reactnity inserted (5) t__.... - . _ . - _ _ _ _ _ _ Fig. 26. Comparison of peak pulsing temperatures from pulses performed at the NSCR and modeled in PARET. The peak fuel temperature is difficult to quantify exactly due to the nature of the experimental data acquired. The experimental temperature was taken from a thermocouple which is located mid-way between the fuel surface and fuel centerline of the instrumented element. Since the measured temperature falls between the predicted centerline and fuel
- surface temperatures for pulses of less than 51.65, the comparison is accepted a3
! reasonable. In addition, the calculated peak fuel temperature is considerably higher than the experimentally measured temperature and therefore is conservative from a safety standpoint. To relate the tnermocouple temperature to the predicted peak temperature in the 1 element, a curve fit was used to develop a functional relationship for the ratio of the peak j fuel temperature to peak thermocouple temperatute as a function of peak centerline temperature. A plot of this for the data from Figure 26 is shown in Figure 27. ; l
46
, 1 65 1 5h 33 te l
EE i 1 55 ff
- 2 jy 1.50 145 j
U
$ .g 140 E5 135 300 00 350.00 400.00 - 450.00 500 00 550.00 600.00 Peak Fuel Centerline Temperature ('C)
Fig. 27. Ratio of the peak fuel temperature calculated for the instrumented fuel rod using
- PARET to the peak thermocouple temperature measured for pulses performed at l the NSCR.
When a second order polynomial in temperature is fit to this data the result is 2 R(Tct)=0.4503 + 4.149x10-3 Tct -3.6857x10~6 T ct (7) where: R(Tet) is the ratio of the peak centerline temperature to the peak thermocouple ; l temperature ; j Tct is the peak centerline temperature of the instrumented element as calculated i by PARET. j l 1 1 A plot of the power and energy obtained from PARET for a sample pulse of $1.65 is i shown in Figure 28. The temperatures calculated from this pulse are shown in Figure 29. These show the power, energy, centerline temperature, fuel surface temperature, and l cladding temperature for the instrumented element for each pulse. 1 i t f { i i i I
47 N h o 4< i i i
, 4 g n ,D 8
'O c) S
~
I'
, d .5 i c
' o i g 3 l' d
~
u & m
$E
- c. u]
i e a g l l 6 t O i e O T l , W
' 4C y t
l u? S fl
, oeu g l
I b F m q C l
~
d d I i e l 8 G 9 wh
-==- s o i
g D @ N O d C l U
'O C
,- c3
~
~
d y U t O
, , i 9 4 o .
9 9 9 9 9 M ( o o o o o N o o o o . t ca (o v N
.b0 (cas-gly) XHJaug *(ggs)Jamod l
i
. . - ~ ... . . _ . . . . . . . . . . - . . . . - - - - - - - - - . - - - - - - - - - - -
500 . . - - - > > ' ' ' ' ' 400 - - h, r Centerfir.e G 300 -i
- - - Fuel Surface -
E 3.
\s
- - Clad e 's E \ s m -
8. e
's s i
1 g 200 -
's. '
r s_
. ~~__
l ~~ ' ' -_. __ 100 - l - r _______
' ' '~ ' ' ' ' ' ' ' ' '
0 ' ' t 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 Time (seconds) Fig. 29. Fuel centerline, fuel surface and clad temperatures for the instrumented channel from a $1.65 pulse at 1400 ImVD as calculated by PARET. . oc o n
. . . , _ _ . . . - . _ . . _ . _ _ _ - . . _ . . . . - . . . . - . - . - - - - ~ - _ . . . . - - . . . . . - - . . _ . - - . . . _ . . . . . . . . . . . . _ . - . . _ . . . . . . _ . _ _ _ _ . _ , - - . . . _ . . . _ . . . . _ . _ . . . . . _ . _ . . _ . . . _ . . . . _ _ _ . . , _ _ . _ _ . . =
i l 49 Core Average Temperature Since the reactivity of TRIGA fuel is extremely sensitive to temperature, the core average temperature must be calculated to determine the appropriate neutron cross sections to be used in the diffusion calculations. The procedure for calculating the core average temperature begins by estimating the core average temperature at the desired power level and then running DIF3D to find the critical control rod positions. Next, the axial power distribution for an " average" rod from the DIF3D case is normalized to an average power l of 1.0 and input to NCTRIGA which is then run at the desired power level to determine the temperature distribution within the rod. if this rod has a radial peaking factor that is 1.0, its temperature distribution is representative of the average fuel temperature distribution i for all of the rods in the core. To determine the core average fuel temperature from the temperature distribution of the " average" rod, the average of the axially varying centerline temperatures and the average of the axially varying fuel surface temperatures are used in the following formula which is used in NCTRIGA to calculate the average fuel temperature. Tove = (Li + 2*Ta)/3 (8) where T,ve = average fuel temperature i Tun = cverage fuel surface temperature i Te = average fuel centerline temperature 1 If the temperature calculated from the above procedure varies greatly from the estimated initial core average temperature, DIF3D is rerun using cross sections generated at the new temperature to fimd the critical control rod positions, and the above procedure is repeated until consistent temperatures in DIF3D and NCTIUGA are obtained.
I l i l l 50 l CilAPTER IV I l SAFETY EVALUATION I i l The calculation of the normal operating parameters, the determination of the Limiting l Safety System Settings (LSSS) and the evaluation of standard reactor accidents are ! presented in this chapter. The values presented here were generated using the computer {; codes described in Chapter III. Wherever possible, the accepted results from the SAR for the currently licensed power level of 1.0 MW are presented along with the newly calculated results for the 1.49 MW power level. I Reactor Analysis The analyses of the physics parameters that would be affected by the power upgrade are presented in this section. The most important values calculated were the peak fuel temperature and the maximum power density expected during normal operation. It was demonstrated that the power upgrade will not cause the temperature limit on the fuel to be ' exceeded. Also, the maximum power density remains, by a substantial margin, below the level at which departure from nucleate boiling could occur. Other values that would need to be included in a safety analysis report amendment for the power upgrade are also included in this section. These include the control rod worths and the shutdown margin. The first step in analyzing the core with the upgraded power 1.49 MW was to determine the core average temperature. It was found to be approximately 550 K. This is somewhat higher than the core average temperature of 500 K found for a power level of 1.0 MW. The increase in the core average temperature decreases the reactivity of the fuel by about $0.50. Therefore, the critical control rod positions at 1.49 MW will be more withdrawn than those at 1.0 MW. Axially averaged power distributions at 1.49 MW and 1.0 MW are shown in Figures 30 and 31 respectively. l
! 52 t l l I , Relative Power 1.72357 1.60867 1.49376 1.37886 1.26395 1.14905 1.03414 0.91924 j 0.804335 0.68943 l 0.474525 0.45962 O.344715 0.22981 l l 0.I14905 i Fig. 31. Relative power distribution of the NSCR at 1.0 MW calculated with DIF3D. l i 4 1 1 i
53 It can be seen from the two previous figures that the power shape will not change significantly in upgrading the power level of the NSCR. The peak to average power at 1.49 MW is slightly higher than that at 1.0 MW, but the overall distribution is quite similar. The highest steady state power density calculated for the NSCR at 1.49 MW was 26.82 kW/ element. Since the Torrey Pines TRIGA Mark 111 has been licensed with a maximum power density of 32 kW/ element, and the NSCR pool is six feet deeper than that at Torrey Pines, the maximum power per element for the upgraded power level is acceptable.2 The highest peak to average value calculated for the rod powers at 1.49 MW was 1.62. To reach the 32 kW/ element limit, the peak to average value would have to reach 1.93. This leaves a large margin as long as the core configuration is not significantly changed. When NCTRIGA was run for the hottest rod in the core at 1.49 MW, only nucleate boiling was calculated to occur. Steady state calculations in PAIWT produced a minimum burn-out ratio in excess of 4 at any axial location along the hottest channel. Although this leaves less margin than the minimum burnout ratio of 6 calculated for 1.0 MW, it still _ provides for extremely safe conditions. The axial power distributions for the hottest rod in the core at 1.0 MW and 1.49 MW are shown in Figure 32. It can be seen from this figure that the axial distribution will not change significantly with the power upgrade. Since the control rods are withdrawn further at 1.49 MW than at 1.0 MW, the peak power in the rod occurs at a slightly higher position. The values of the peaking factors and overall shape of the axial distributions changes very little with the increase in power. The peak fuel temperature calculated at 1.49 MW is 593 C and is 150 C higher than the peak fuel temperature calculated at 1.0 MW. The coolant mass flow rate in the hot channel increased to 284 kg/m2 s at 1.49 MW from 244 kg/m2s at 1.0 MW. l l l
54 Figure 33 depicts temperature values calculated from NCTRIGA for the " average" channel in the core for both 1.49 MW and 1.0 MW. The values shown there are i i representative of core average values that are expected during operation at each of the two , i l power levels. t 2.5 ' ! 2 ' l .- . i e j j1.5 ,,,, l 5 I I 1 H ! 0.5 - - - - - - 1.0 MW , 1.49 MW j i 0 i 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 i Atlal licight (m) Fig. 32. Axial distributions for the hottest elements at 1.0 MW and 1.49 MW calculated with DIF3D. I I l [ . l i
55 I
' ~
450 --+-l fuel cl (139 MW) ;
, +Twall (1.49 MW) I 400 e Twater (1.49 MW)
-M-Truel cl (1.0 MW) >{'
350 TwaH (1.0 MW) 4 o v
I --*--Twater (1.0 MW) a C 300 m
e i250 2 E200 E
@ 150
- : : : : : : : : : A = = = =
100 i i 50 e e e 0 ?---*---
- I" #
m
, e e I
0
' I' O 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 Asiallleight (m)
Fig. 33. Fuel centerline, wall, and water temperatures for the " average" channel calculated with NCTRIGA at 1.0 MW and 1.49 MW. It can be seen from Figure 33 that the fuel centerline temperature will increase by approximately 100 C, but the fuel surface temperature, clad temperature, and water temperature will hardly increase. The centerline temperature increase is due to the increase in the power generated per element. The lack of a temperature increase at the fuel surface and in the coolant is caused by the increase in the mass How rate. The mass flow rate at 1.0 MW was calculated to be 207 kg/m2s and the mass flow rate at 1.49 MW was calculated to be 240 kg/m 2s for the average channel.
- The calculated heat Duxes for the " average" channel are shown in Figure 34 and the l
calculated heat transfer coefficients and Reynolds's numbers are shown in Figure 35 for the " average" channel at both 1.0 MW and 1.49 MW. Again, these parameters increased due to the higher power level, as expected. l t
56 l ! 600000 500000 ^ .
-+-Heat Flux (1.49 MW)
+ Heat Flux (1.0 MW)
.~.
$ 400000 E ,.
E 300000 C j 200000 ( 100000 l 0 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 Axial lleight (m) l Fig. 34. Heat fluxes for the " average" chanr:1 calculated with NCTRIGA at 1.0 MW and 1.49 MW. i 12000 -
-+--H T Coef (1.49 MW) ,
y -E-Reynolds (1.49 MW) ' E 10000 H T Coef (1.0 MW) I g X Reynolds (1.0 MW)
;g 8000
'E ! .
I
) 6000 U3 1
$E PA-- -h--A I gA d %..
g 2000 lE O O 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 Axial lleight (m) L l Fig. 35. Heat transfer coefficients and Reynolds's numbers for the " average" channel calculated with NCTRIGA at 1.0 MW and 1.49 MW. I
57 All of the above parameters were calculated to change as expected with the power upgrade. Therefore, the NCTRIGA model is believed to perform satisfactorily for the NSCR. No cooling problems should be encountered while operating at 1.49 MW. Control rod worths are of great importance for maintaining reactor control. Table 7 l compares the calculated control rod worths at the core average temperatures corresponding l l to 1.49 MW,1.0 MW, and 300 W, along with measured values from January 1995. The rods are worth slightly less at 1.49 MW than at 1.0 MW, but the overall worth is only l $0.25 less. TABLE 7 ! Control Rod Worths Measured at the NSCR and Calculated with DIF3D Calculated with Calculated with Calet'Uted with Measured Temperature Temperature Temperature 1/6/95 l Corresponding to Corresponding to Correspopding to ; i 1.49 M W l.0 M W 300 W l Configuration keg dk Worth ky e dk Worth kg e dk Worth Worth (S) ($) (S) ($) l All Rods In 0.928 0.931 0.941 l SSI Out 0.946 0.018 2.63 0.950 0.019 2.69 0.959 0.019 2.68 2.84 l SS2 Out 0.941 0.012 1.78 0.944 0.013 1.82 0.953 0.013 1.82 1,91 SS3 Out 0.948 0.020 2.81 0.952 0.020 2.87 0.961 0.020 2.86 2.89 L SS4 Out 0.961 0.033 4.73 0.965 0.034 4.82 0.974 0.033 4.78 4.75 RR Out 0.938 0.010 1.40 0.941 0.010 1.40 0.950 0.010 1.41 0.82 TR Out 0.947 0.019 2.68 0.950 0.019 2.67 0.959 0.019 2.67 3.24 l Total 16.02 16.28 16.21 16.45 l I Another required value calculated for the NSCR was the shutdown margin. The first , computation was to determine the excess reactivity at 300 watts. To calculate this, critical control rod positions were found with the transient rod fully withdrawn and a core average , temperature of 300 K. Next, the keg was calculated with one of the rods fully withdrawn, 1
58 and the remaining rods still in their previous positions from the critical configuration. The change in kerr for each case represents the remaining worth of the rod ofinterest. Table 8 shows the calculated worth remaining for each control rod when they are in a critical configuration with a core average temperature of 300 K. TABLE 8 Excess Reactivity at 300 Watts Measured at the NSCR and Calculated Using DIF3D Remaining Worth ($) Calculated Measured 1/6/95 SSI Out 0.9639714 1.16 SS2 Out 0.642076 0.85 SS3 Out 0.8706021 1.32 SS4 Out 1.4182713 2.24 RR Out 0.2628406 0.52 Total 4.1577614 6.09 The difference in the calculated and measured values could be the result of choosing a different critical configuration for the calculation than was used in the measurement. The details of the experimental critical configuration were not available. The calculation of the shutdown margin is defined as the total control rod worth less the sum of the following: the core excess, the highest worth non-secured experiment, the l most reactive control rod worth, and the worth of any non-scrammable rods. The i calculation of the shutdown margin for a core average temperature of 300 K is shown m i i l Table 9. The temperature of 300 K was chosen as the ambient temperature for all cases 1 i 1 run. The shutdown margin is calculated at this temperature since it is the coldest and most J reactive the core will ever be. The shutdown margin here is smaller than it would be at any 1 other temperature. I 1 l
1 1 ! l ! l i 59 I l l Since the shutdown margin is calculated at ambient temperature regardless of the maximum power level, the shutdown margin will not change with the power upgrade. It is presented here because it is a necessary calculation for the SAR. amendment and it has been calculated using the new methods. l l TABLE 9 l Shutdown Margin Measured at the NSCR and Calculated Using DIF3D ! l Calculated Measured l 1/6/95 ! fotal Control Rod Worth 16.21 16.45 Core Excess at 300 Watts -4.16 -6.09 l liighest Worth Non-Secured Experiment -1.00 -1.00 Most Reactive Rod Worth -4.78 -4.75 Regulating Rod Worth 2.67 -0.82 Shutdown Margin 3.60 3.79 l l Evaluation of the Limiting Safety System Setting ) The Limiting Safety System Setting (LSSS) is a temperature that will initiate a scram when measured by the thermocouple in the instrumented fuel element. The LSSS is established to ensure that the safety limits will not be surpassed under any condition of operation. Although the LSSS is not expected to be significantly affected by the power i upgrade, it must be calculated for the SAR amendment. A peak fuel temperature of 950 C j has been established for all modes of operation with FLIP fuel and it provides a minimum safety margin of 200 C. Since the fuel temperature can only be measured in the .
~
instrumented fuel element, the difference in the measured temperature and the peak fuel temperature in the core must be conservatively calculated. The LSSS is established in the l c SAR as the lowest temperature at the thermocouple location when the peak fuel
l l l 60 i temperature in the core is 950 C.4 Analyses were performed to establish the LSSS for both steady state operation and reactor pulsing. The Safety Analysis Report allows the instrumented fuel element to be located at any one of eight locations in the core that are adjacent to the central bundle, excluding the corner positions. These locations are shown graphical y in Figure 36. The position of the IF is allowed in these eight locations so that fuel shuffling can be performed. As a result of this procedure, all eight locations must be examined, and the location that will produce the lowest temperature in the IF rod must be used in determining the LSSS. This typically results in the LSSS being set low, but this policy must be followed in case the IF rod is ever loaded in the position with the lowest power. Since the power distribution is different for steady state and pulsing conditions, these two cases are considered separately. LSSS for Steady State Operation DIF3D was run for a power level of I A9 MW using the fuel burnup from the reactor records as of January 1995. The core average temperature at 1.49 MW was calculated to be 550 K. NCTRIGA was then run using the axial power distribution from the rod with the highest radial peaking factor. The power level was increased in the NCTRIGA model until a peak fuel temperature of 950 C was reached. This power level was 2.6 MW, and the corresponding average core temperature was found to be 661 K. DIF3D then run using cross sections generated at 650 K, and the axial power distribution for the hottest rod from this case was input to NCTRIGA. The core power to produce 950 C in this case was 2.5 l MW, and the core average temperature was found to be 647 K. Thus, the DIF3D job was concluded to be run at the correct temperature, and this case was used to establish the LSSS for steady state operation. A plot of the peak fuel temperature as a function of reactor power, as calculated with NCTRIGA, is shown in Figure 37. f i \
61 A B C D E F 4 o n M tk& o i MiRBRM o c 8383g: 4 , . o 1. q: jp es? .
^
e.. q i,1 $ $ fijk oOpeL*JEg 1 - f .i ; - L.. .t . t i I Posrible IF Rod Locations FLIP Fuel [, Graphite Redector Fig. 36. NSCR core map with possible IF rod locations. 1
-62
, 966 ; l 964 , G l y962-5960 3 !
$958
- m. ;
l ! 556 - ! t l w X954 j952 n 950 - , d ' ' ' 948 2.5 2.51 2.52 2.53 2.54 2.55 Core Power (MW) , Fig. 37. Peak fuel temperature as a function of reactor power calculated with NCTRIGA. l l The peak temperature in the lowest power IF position was found to be 812*C at a core ! power of 2.5 MW. Since the highest observed ambient temperature for the pool is 37 C, the peak temperature rise of this instrumented fuel element would be: 1 . l ATT c = 812 C - 37 C = 775 C (9) 1 1 In steady state mode the radial temperature distribution through a fuel element is not ( highly sensitive to the within rod power density distribution. Thus, the radial temperature t l j distribution which was calculated for the Puerto Rico Nuclear Center reactor is applicable ! to the NSCR. For this case the ratio of the maximum temperature, which occurs at the 4 centerline of the fuel element, to the temperature at the thermocouple location was cited in the NSCR SAR to be 1.16.4 Thus, for the element discussed above, the thermocouple I should read i
e.%a f ;Aa e4~ Z. 4E-ha--44 J
- Ml.+,., 4 J. .--M.J= 4.4 in 44m& 4 *+ae rh _+65.n..d 4._ 5.. 4A _m, _ . 2 ,_A,_a_ 4 4 , 3 63 775 C
=668 C above ambient (10) 1.16 Thus, for the FLIP core, the LSSS is LSSS = 668 C + 37 C = 705 C (11) i This is lower than the value of 728 C presented in the previous SAR. This value is not dependent on the power level that the core operates at. The LSSS for steady state !
L operation is based on overpower conditions, to assure that fuel damage will not occur. The LSSS and Maximum Allowable Reactivity Insertion for Pulsing Any scram that occurs due to temperature during a pulse needs to prevent the safety limit from being exceeded by compensating for the fact that the temperature of the fuel ; continues to increase for several seconds after a pulse. Control over the maximum temperature incurred during pulsing is achieved by limiting the amount of reactivity that can be inserted during pulsing so that the limiting temperature is not reached. It is, however, accepted procedure to set the LSSS as the thermocouple temperature to assure peak core temperature no greater than 950 C during a pulse. The minimum LSSS for pulsing was determined for the thermocouple location adjacent to the central fuel bundle that establishes the largest value for the power ratio (P/PTc). This power ratio is the ratio of the maximum core power, P , to the power at the thermocouple l 1 i location, Pic. To find the reactivity insertion that will produce a peak temperature of 950 C, PARET was run using the power distributions for the hot channel generated by DIF3D based on l prepulsing conditions. Since it is well known that the end of core life pulses produce the highest temperatures for a given reactivity insertion, only cases at end of core life were i I
l l l 64 examined. The radial power distribution from the prepulsing I)IF31) case is shown in Figure 38. Various pulses were modeled until a peak core temperature of 950 C was obtained. The reactivity insertion required to produce a peak core temperature of 950 C was found to be $2.66. This compares well with the previously accepted value of $2.48.4 A plot of the pulse power and energy is shown in Figure 39, and a plot of three temperatures from the hot channel is shown in Figure 40. PARET was then run to determine the temperature response at the lowest power instrumented location. This rod had a peak temperature of 800 C during the $2.66 pulse. The temperature of the thermocouple was then determined from this. Evaluating the numerical fit from equation (7) of Chapter 111 for a peak temperature of 800 C to find the ratio of the peak to thermocouple temperature rise gives, 4 R(800 C)=0.4503 + 4.1419x10-3 x (800)-3.6857x10 x (800)2 (12) R (800 C) = 141 (13) This gives a thermocouple temperature of 800 C-37 C
+ 37 C = 541 C (14) 1.41 which is taken to be the LSSS for pulsing operation. This is very close to the value accepted from the previous SAR of 544 C for an all FLIP core. A list of temperatures and l peak powers for several pulses near the pulsing limit are shown in Table 10. It can be seen i from this table that the fuel temperature and peak reactor power are quite sensitive to changes in the reactivity inserted. A plot of the calculated peak fuel temperature and thermocouple temperature is shown in Figure 41. This shows that the peak fuel temperature approximates a linear function of reactivity insertion.
1
( . l i 65 l l l i 1 i 4 l Relative Power 1.69567 l 1.58262 5 1.46958 1 1.35653
, l.24349 l
1.13044 1.01740 0.904355 0.791311 h ,. f 0.678266 4 N ~ 0.565222 0.452177 E u, 0.339133 w;; i 0.226089 i I i 0.113044 i I l l Fig. 38. Relative power distribution calculated by DIF3D for the prepulsing case. 1 i 1 1 4
#*#*'4 d-h-eEimu %4,., g 4 4_
i
-f 66 i
M I y 6 l a 8 i , > O. g e i T"" sp) i "D
' 4 l R d i i O O
' W f
~
tE 80 O k N i M' i SE :
, w
, 1 g 9 0
~ l & .
g I O i 4 l f i I b 9 0 -; I i ' s O h ! o ? D
- R' o ;
i
- > O o i m.
O O ~ C3 o
- E o
l l t i v. O P o3 c5 e 1 h t l l 0 o b n = y s f. I Q $ 6 t 1 4 .
-g f I
'M I
f U g 4 O
- e. %
O . Ch M t ' ' e
- o. - .%
O o O o, g o 4-d 6 e O O o O O c O o 6 O O O l D 4 o m - (335-/AW) ASJaug *(g}q)ja*0d i I i 1
'N -+_____m
i
?
i t 1000.0 . . . i i i i i . . . .. . .
}
t 800.0 -l - l
,'i .
, t i ,
i i
' Centerline !
G 600.0 -
'i - - - Fuel Surface
[ U ' - - aad ! 3 ' i t' 2 i il- ', '
- 8. 's s ,
e
# 400.0 -
-N -
, s :
's f i
200.0 -
~~_ -
~
s ___ ~~ __ ____ i s - -__
~
0.0 , 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 Time (seconds)
. i Fig.40. Fuel centerline, fuel surface and clad temperatures calculated by PARET for $2.66 EOL pulse. ,
" l 4
9 i
~(
68 l l TABLE 10 ' Temperatures and Powers Calculated Using PARET for Pulse Insertions at EOL Reactivity Insertion Power (MW) Peak Temperature Thermocouple ($) ( C) Temperature ( C) 2.00 1474 663 321 2.25 2365 770 384 2.50 3504 878 469 2.66 4372 949 543 3.00 6668 1103 813 1200
>l 1000 i
800
,d l 4 -
T < l
$ 600 ' '
c-E "V-g, O Peak Fuel
--E--Thermocouple 200 0
2 2.2 2.4 2.6 2.8 3 Reactivity Insertion (S) Fig. 41. Peak fuel and thermocouple temperatures calculated with PARET for prompt reactivity insertions at EOL. , f
69 The limiting values for pulsing are shown in Table 1i for the newly calculated values and the previously accepted values. Again, these values are not dependent on the normal operational power level of the core. The difference between the values presented in the SAR and the valued calculated here is the result of the new methods. TABLE 11 Limiting Safety System Setting Values for Pulsing I Maximum Allowed Pulsing LSSS Reactivity Insertion From June 1979 SAR 544 C $2.48 Newly Calculated 541"C $2.51 Since the LSSS for pulsing is lower than the LSSS for steady state operation, it will be used as the LSSS for all modes of operation. To ensure that this limit would not be exceeded during normal operation, the highest possible thermocouple temperature that can i be reached during steady state operation was investigated. The highest power IF rod I location gives a peak temperature of 553 C for a core power of 1.49 MW. Using the same ratio of the centerline temperature to the thermocouple temperature used earlier gives a thermocouple temperature of: 1 l 553 C- 37 C + 37 C = 482 C (15) 116 The margin of 59'C between the peak thermocouple temperature and the LSSS should be ! adequate for normal operation. It was calculated that a core power of 1.83 MW would have to be reached before the thermocouple temperature would reach the LSSS of 541 C. I
l L 70 i i Accidra Analysis l The four accidents analyzed in the current SAR for 1.0 MW are analyzed here for ; 1.49 MW. These accidents are: accidental pulse at full power, control rod nm-out, loss of I coolant accident, and design basis accident. The accidental pulse at full power is actually ! I an extremely unlikely accident, since the pulsing mechanism is locked out when the ! reactor is at power. The objective of the analysis of the accidental pulse at full power was to find the : l reactivity insertion required to achieve a peak fuel temperature of 950 C when the reactor l i is operating at 1.49 MW. Since it was known that end-of-life pulsing would produce the l t highest fuel temperatures, only pulses at end-of-life were analyzed. PARET was run using l a power distribution for the hot rod from a DIF3D job run with the transient rod inserted and a core average burnup of 3000 MWD, which is the stated end of core life.d The , i prompt negative temperature coefficients for 3000 MWD burnup were also used. The ! reactivity insertion that produced a peak fuel temperature of 950 C was $1.96 which is much lower than the SAR value of $2.42 for an insertion from 1.0 MW. The power and I t energy produced during this pulse are shown graphically in Figure 42. The temperatures i produced during this pulse are shown in Figure 43. l Although this accident is not likely due to safety lockouts, the analysis is required for licensing. It is possible, however, that an experiment with a high negative worth could detach itself from the core during operation. This would have a similar affect to ejecting l the transient rod at power. Since $1.96 could be inserted with the transient rod before temperature limits on the fuel are exceeded, the $1.00 limit on the worth of a non-secured experiment will still allow for safe operation. The control rod run out accident is a another transient that must be analyzed. It involves driving out all four shim safety control rods in gang at their maximum speed. This accident is analyzed at the end-of-life because the reactivity feedback is smallest at 1 _ _ , - ._ _ -- _ _ ~ , - - , - . -. m
v f t 2500.0 . . . . . . . , , i t f
.2000.0 - -
i Power i g . . . . . - Energy , i i 1500.0 - - i ll2 u ' di I E 1000.0 - t . 500.0 - 0.0 ' '-- ' ' ' 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 r time (seconds) Fig.42. Power and energy calculated by PARET for an accidental pulse of $1.96 from full power at EOL. y .-- - - - - . _ . - - - - - . - - . - - - - . - - . - - - - . - . . . . - - - - . . - - - - - - _ . - - . - - - - - - - - - --.--.---=------r. ..2 , - . . w. s- -u~.+ . - > ----e,-, . ,.-n,-- . . e ,,,a-r r,n e ,m,m--,-* m w
i 1000.0 , , , . . i i i i . . . i i 800.0 - - i i l 4 } i 4 i G 600.0 -
'i centertine -
y F
', - - - Fuel Surface
{a . i - - Clad i s , E i 's, g, i g a % g 400.0 -l 'N - t I s t 's ~,
.i
,i s ,
i, _ t
, ~
200.0 t - - t,r ~ _ x ____
~
_ _ _ _ mme 0.0 ' ' ' ' ' O 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 , Time (seconds) Fig.43. Fuel centerline, fuel surface and clad temperatures calculated by ; PARET for an accidental pulse of $1.96 from full power at EOL. y l M :
73 l that time. Since it was not entirely clear from which power level the accident should be analyzed, two analyses were performed. The first case investigated running out the control rods from an initial core power level of I Watt. To analyze this accident, a cold critical configuration was found with a core average burnup of 3000 MWD. The integral worth of the shim safety bank was also determined using D'J3D and is shown in Figure 44. The power distribution for the hottest rod in the core and the bank rod worth were input to PARET. The overpower trip setpoint was fixed 1.86 MW which corresponds to 125% power. The maximum rate at which the current control rod drives operate is 0.19 cm/sec, but new rod drives that could have a withdrawal rate as high as 0.59 cm/see will soon be installed. Therefore, this case was analyzed as a worse case scenario. The maximum fuel temperature produced in the transient was 99 C, which creates no safety concern. A plot of the power and energy from this transient is shown in Figure 45, and the temperatures are shown in Figure 46. The power shape looks much like the shape encountered from a pulse insertion. This was caused by the reactor being placed in a prompt supercritical configuration with over $1.00 in control md worth inserted before the scram occurred. The second control rod run out accir'.ent analysis was performed from an initial power level of 1,49 MW. For this case, the peak core temperature was 516 C, which is well below the temperature limit on the fuel. A plot of the calculated power for this case is shown in Figure 47, and a plot of the calculated temperatures is shown in Figure 48. The power shape from this accident looks considerably different than the shape from the previous accident. This is because the reactor is in a supercritical condition when the rods are withdrawn from 1.49 MW, as opposed to the prompt supercritical condition that is eventually reached when starting from 1.0 Watt. When starting from 1.49 MW, the scram occurs in less than one second from the beginning of the transient as opposed to the five second time interval experienced when starting from 1.0 Watt. It should be noted that the power produced in the transient actually exceeds the trip setpoint of 1.86 MW. This is
r i 74 t l 12 , i l l ? 10 ? i l t t B l> l ! i R
- v. ;
.c >
r j 6 l i .x
- c l .!!
i i - I 4 ; t 2 . i J 0 l om 10m 20.00 mm 4 00 50.00 som 70.00 som 90 2 100m ) Percent Withdrawn i i Fig 44. Integral rod worth for the shim safety bank calculated with DIF3D. i i 1 l 1 i t j
4 4 25 . . . . . . . . . . . . . . 20 - c-o fg 15 -
..-- Energy Power -
E u e E 10 - t 5 - - 0 ' ' ' O 1 2 3 4 5 6 7 8 9 10 11 12 13 -14 15 i Time (seconds) Fig. 45. Power and energy calculated by PARET for a control rod run out from 1 Watt at EOL. u
*n ?~ w- ------ _-w a*+ wee-ee.* e--eow h e-e m-g- g-- M *'=***->+e==+av w -m e e- *we*> =uws-*'N*=rtwm-+e' w-e.e wrer a-- e um vr- W"wuh=fraiww+e-'wewe==-*2"+'-- m'* -*Fa' *----w* - - - - - - - - - - - - - =-- - - - - - - - - - - - - - "
100.0 . . . . - - i . i ' > > - ' I s_ s 80.0 - fs s f 's , i 's s , [ 's s - G I 's s '
- P
-o 1 's s
; i ~s s' ~_
g 60.0 - l ~~ _- j
~.
8 l 5 e i l l I I Centerline
- 40.0 -
ji - - - Fuet surface l - - Clad i 20.0 O 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 Time (seconds) Fig.46. Fuel centerline, fuel surface and clad temperatures for the hot channel calculated ' by PARET for a control rod run out from 1 Watt at EOL. s e
. . . .--_.-,,.___.,.--.,.__m. _ ______..__,.,_.,-.m.--, - . - . _ _ , , . . . . _ . . . . . - . . . . . _ _ _ - - . e,,4 ....m_--
1 2.0 , , . . . . . . . . . . . . . . . . . 1.5 1 F
- E r 1.0 -
4 2 i 0.5 - - i t 0.0 ' ' ' O 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 Time (seconds) Fig.47. Power calculated by PARET for a control rod run out from 1.49 MW at EOL. O e k 4 1
- _ _ _ _ _ _ _ _ _ _ . m.,e.-_ ---..2-w., -.,,--,i..-,.=-
600.0 . i i i i i i i i i i i i i Centerline
- - - Fuel Surface 500.0 -
- - clad 400.0 -
G F 3 , o aB 300.0
-t i
E .I , e s
# 1s i ',
N 200.0 3 i - 4
---_______________~_~ ---
_ __ _ -~ _
=
100.0 - 0.0 O 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 Time (seconds) Fig.48. Fuel centerline, fuel surface and clad temperatures for the hot channel calculated by PARET for a control rod run out from 1.49 MW at EOL. y oo
r i 79 caused by the 0.2 second delay input to PARET to account for the reaction time of the control rod release mechanism to the trip detection. Although a control rod run out accident from 1.0 MW is not necessary for a license l amendment, it was analyzed for comparison with the control rod run out from 1.49 MW. The peak temperature reached in this accident was 402*C. Again, this temperature is well below the temperature limit on the fuel. A plot of the calculated power from this accident j is shown in Figure 49, and a plot of the calculated temperatures for the hot channel is l shown in Figure 50. l
- The loss of coolant accident (LOCA) is scenario in which the reactor pool is accidentally drained of water through a 10 in. diameter hole at its bottom. Heat from the fuel is then removed only by the natural convection of air. Since the reactor will scram if the pool level becomes too low, only the decay heat must be removed. The limiting value for this scenario is the maximum temperature that FLIP fuel can tolerate in air before the cladding fails, which is 940 C. A thorough analysis of this accident was performed by GA and presented in the SAR. The conclusion of their analysis was that FLIP fuel will not exceed 940 C during a LOCA if the maximum power density does not exceed 28 kW/ element and reactor operation is limited to 70 hours per week. The analysis of the NSCR at 1.49 MW gave a maximum power density of 26.82 kW/ element, and the NSCR no longer operates for more than 70 hours per week. Thus, the power upgrade will not create a LOCA safety hazard.
The design basis accident is defined as the loss ofintegrity of the fuel cladding for one ; l fuel element and the simultaneous loss of pool water. During the accident, fission products l are released from the fuel. Again, this accident was thoroughly investigated in the current SAR using a limiting power density of 28 kW/ element. As stated above, the maximum l power density for 1.49 MW is 26.82 kW/ element, so the previous analysis remains valid. I Doses to persons inside and outside the facility were calculated with and without the pool j t l f
4 i i i 4 2.5 . . . . . . . . . . . . i . . . . . . p t 2.0 - r I r i 1.5
- i. F >
E . t i s n. i 1.0 - I t
- 0.5 - -
i I i 4 0.0 O 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16.17 18 19 20 : Time (seconds) l Fig.49. Power calculated by PARET for control rod run out from 1.0 MW at EOL. : oo t t i t ____.___.________m._.___.m._ _ _ _ _ _ _ _ , _ _ , _ . . , . . _ _ , . _ _ _
,,__.~,,m_..,,_.,_..,_.m_, , _ , _ , , , , , ,,_, ., ..,,,_ , ,, , .., . _ , _ _ _ __,_.m_,,., , _ . _ , , , , , , , , . , , _ , ,, , ,, ,,
500.0 . . . . . . . . . . . . . . Centerline
- - - Fuel Surface
- - c ad 400.0 _ _
G 300.0 - - t' 3 H
& i
-e r g 200.0 4 I
~ _ _ _ - ,
______~~~_ ._ ~~. ______- __ ~~~ 100.0 -
.~___' '
0.0 O 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 Time (seconds) ! Fig. 50. Fuel centerline, fuel surface and clad temperatures for the hot channel calculated l by PARET for a control rod run out from 1.0 MW at EOL. g i l l t
l 82 water remaining as well as with and without the building ventilation system operating. l The results of these calculations are presei.ted below in Table 12. l TABLE 12 Summary of Radiation Exposure Following the Design Basis Accident WBGD Thyroid Dose Ventilation Pool Water System Maximum Exposure to 3.5x10~3 mr -- operating remaining Population Outside Building 1.4x10-2 mr 3.2 mr operating drained 0.84 mr -- shut down remaining 1.75 mr 10.5 r shut down drained Maximum Exposure to 3.6x10-3 mr -- operating remaining Operating Personnel in 2.1x10-2 mr 18 mr operating drained One Hour After Release 1.75 mr -- shut down remaining 4.2 mr 49r shut down drained For a significant exposure to occur, the pool water would have to be lost at the same time the ventilation system fails and a fuel rod ruptures. Even in this case, only personnel who remain inside the facility would be effected. As stated in the SAR, the design basis accident presents no realistic hazard of consequence.
83 l CIIAPTER V , CONCLUSIONS The core parameters presented in the SAR used to license operation at 1.0 MW have ; t been compared to those for operation at 1.49 MW. It was shown that the maximum power density in the core would not exceed 26.82 kW/ element during steady state operation. Since other TRIGA cores have been licensed to operate with power densities as high as 32 kW/ element, and the margin to departure from nucleate boiling was adequate, 27 kW/ element should be acceptable at the NSCR. Other physics parameters investigated such as the relative power distribution, control rod worths, and shutdown margin would not be significantly effected by the power upgrade. Therefore, the analyses presented earlier demonstrate that the NSCR can be safely operated at a power level of 1.49 MW, upgraded from 1.0 MW without additional modification to the core or control rod drives. The peak fuel temperature calculated for s: ady state operation at 1.49 MW was 593 C. This is well below the temperature limit on TRIGA fuel for this facility of 950 C. The actual minimum temperature at which cladding failure could occur is 1150*C. The Limiting Safety System Setting changes little with the new methods. The previous setting corresponded to a thermocouple temperature of 544 C, and the newly calculated value is 541 C. The highest possible thermocouple temperature produced during steady state operation at 1.49 MW was found to be 482 C. This would allow for normal steady state operation at 1.49 MW without the possibility of a scram occurring due to the LSSS, which was fixed by pulsing conditions. The reactivity insertion limit for reactor pulsing from 300 Watts was calculated to be $2.66, which is slightly higher than the licensed value of $2.48. This difference is mostly due to the slightly different rod peaking factor and j temperature coefficients for the fuel from the newly applied methods. In practice, reactivity insertions for pulsing typically do not exceed $2.00, so the pulsing limit does not i restrict operation. l r , e ,
84 l l The accidents analyzed in the SAR were investigated in this thesis. It was found that no greater risk would be incurred by upgrading the reactor power to 1.49 MW than was calculated for 1.0 MW. In the analysis of an accidental pulse at full power, it was found that a reactivity insertion of $1.96 from a core power of 1.49 MW produced the limiting temperature in the fuel. This is lower than the limiting insertion of $2.42 stated in the SAR ' for an accidental pulse from 1.0 MW. This accident should not occur due to safety lock outs. However, a reactivity insertion of $1.00 could occur if an experiment with a high ; negative worth was suddenly detached from the core at full power. Since the highest ! L worth for a non-secured experiment in the core is $1.00, no temperature limits would be exceeded even if the experiment was detached. The analysis of the control rod run out accident from 1 Watt produced a peak core j temperature of only 99 C, and the analysis from 1.49 MW produced a peak core , temperature of only 516 C. Both of these are well below the temperature limit on the fuel , of 950 C. This analysis was performed using the maximum withdrawal rate that could occur even during a malfunction of the control rod drive mechanisms soon to be installed ! in the NSCR. j l ! t l The analyses of the loss of coolant accident and the design basis accident were not ;
- }
i recalculated in this thesis since the maximum power density at 1.49 MW will not exceed ; the maximum power density examined for these accidents in the current SAR. Both of j these accidents were found to have no significant consequences on personnel or facilities. l , The modem computer models used to perform the above analyses have been verified to , l for accuracy with benchmark cases. After review by an independent party, the upgrade of i f the power level at the NSCR should be allowed based on the calculations of this thesis. ! t
)
i l k
, -- i+- .~._.-. ,, - - - - , - - - -
85 Future Work l Further work to be done in this area includes enhancement of the DIF3D model by adding a more accurate control rod model. This work is currently being investigated by j Mark Bigler.33 It is believed that this model will change the control rod worths calculated in this thesis, making them more consistent with the measured values. However, the effects of these changes on the temperatures calculated in the fuel should be minimal since the control rod positions for the critical cases should not change significantly. 1 l The shutdown margin and excess reactivity calculations should be investigated further when the new control rod model is implemented. The worths of the transient rod and regulating rod did not match well using the current model and are believed to have caused ; the mismatches between the measured and calculated values for the shutdown margin and excess reactivity. Although the NCTRIGA model is believed to produce accurate results, its mass flow rate and temperature calculations have not been thoroughly benchmarked with l experimental measurements at the NSCR. This is due to a lack of flow rate information l from the NSCR. Also, fuel temperatures can only be measured at the thermocouple location which is not at the highest temperature location of the fuel element. More l confidence could be placed in the calculations if the peak temperature in a fuel element - could be measured. i I
86 REFERENCES
- 1. Parish, Theodore A.," Nuclear Safety Analyses and Core Design Calculations to Convert the Texas A&M University Nuclear Science Center Reactor to Low Enriched Uranium Fuel," Nuclear Engineering Department, Texas A&M University, College Station, TX (1995).
- 2. Safety Analysis Report for the Torrey Pines TRIGA Mark III Reactor, Gulf General Atomic Project 1252.0000, San Diego, CA (1970).
- 3. Research Reactor Core Conversion from the Use of Highly Enriched Uranium to the Use of Low Enriched Uranium Fuels Guidebook, IAEA-TECDOC-233, Physics Section, International Atomic Energy Agency, Vienna,(1980).
4 . Safety Analysis Report for the Nuclear Science Center Reactor, Nuclear Science Center, Texas A&M University, College Station, TX (June 1979).
- 5. TRIGA Mark Ill Reactor Hazards Analysis, Gulf General Atomic Project 3886 (Rev. A), San Diego, CA (Feb.1965).
- 6. Safety Evaluation Report on High-Uranium Content, Low-Enriched Uranium-Zirconium flydride Fuels for TRIGA Reactors, NUREG-1282, U.S. Nuclear Regulatory Commission (1987).
- 7. C. I. Costescu, W. L. Woodruff and D. G. Cacuci, "WIMSD4m --- WIMS with Isotopes Capability Basic Information for Users," University ofIllinois at Urbana-Champaign, Nuclear Engineering Department (1992).
- 8. S. I. Bhuiyan, Anisur Rashid Khan, M. M. Sarker, M. Rahman, Z. Gulshan, et al.
" Generation of a Library for Reactor Calculations and Some Applications in Core and Safety Parameter Studies of the 3-MW TRIGA MARK-II Research Reactor," Nuclear ; Technology, March 1992, Vol. 97, p. 253-63. l
- 9. K. L. Derstine,"DIF3D : A Code to Solve One , Two- and Three-Dimensional Finite-Difference Diffusion Theory Problems", ANL-82-64, Argonne National Laboratory, Argonne, IL (1984).
- 10. Smith, R.S.," Comparison of NCTRIGA Results to GA Data, NCTRIGA Input 1 Format, and Revision of NATCON," Intra-Laboratory Memo, Argonne National I Laboratory, Argonne,IL (Sept. 18,1992).
1
87
- 11. C.F. Obenchain,"PARET--A Program for the Analysis of Reactor Transients." IDO-17282, Argonne National Laboratory, Argonne,IL (January 1969).
- 12. Schlicht, Roger, FAX correspondence, General Atomic, San Diego, CA (October,1993).
- 13. Bigler, Mark A., Thesis work in progress, Department of Nuclear Engineering, Texas A&M University, College Station, TX (1995).
l I I
88 ! APPENDIX A l WIMSD4/M MODELING The neutron cross section libraries were developed prior to the writing of this thesis. This work was done primarily by Chien-Hsaing Chen at the Nuclear Science Center. A general overview of the methods used is presented here since the accuracy of this work is crucial to this thesis. The code selected for the generation of few group neutron cross sections representative of various subregions in the con and its surroundings is WIMSD4/m, a one-dimensional neutron transport code capable f rolving for flux distributions using many fine energy groups. WIMSD4/m is an advance, version of WIMSD4. The associated cross section data library has three new isotopes: hydrogen in ZrH, '"Er, and *Er. WIMSD4/m reads 69 i group cross section data for each material present, and creates a spatially and energetically averaged set of cross sections in the seven standard energy groups used in TRIGA analyses for each subregion of the core. These cross sections are stored in the ISOTXS format for I 1 use by DIF3D. Since neutron cross sections are highly temperature dependent, WIMSD4/m must be run for many cases covering the entire operational temperature range for each type of l l subregion in the core for both normal and accident conditions. The composition of the core I l components and techniques used for the WIMS modeling are described below. ! l l Fuel Rod Model i Axially, the fuel element is modeled as 38.10 cm (15 in.) of fuel material, capped on both ends by 8.89 cm each of graphite reflector. The bottom fitting is 8.89 cm of stainless steel; the top fitting is 11.43 cm of stainless steel and aluminum.
l 89 In the x-y plane, each fuel cell is modeled as a square with sides of 3.9520 cm. At the center of the fuel meat is a 0.4572 cm diameter zirconium rod. This rod is surrounded by fuel-moderator material with an outer diameter of 3.48234 cm. The meat is surrounded by a stainless steel clad with a thickness of 0.0508 cm. The rest of the 3.9520 cm x 3.9520 cm cell is filled with water. The atom density for each isotope in each material is listed in Table 13. The end reflector cell has the same outside diameter and thickness for the clad as the fuel cell, but is composed of graphite in a stainless steel clad. A diagram of the fuel rod model showing both axial and radial geometries is shown in Figure 51. Control Rod Model There are three types of control rods in the NSCR: shim safety rods, a transient rod and a regulating rod. To model the control rods in WIMSD4/m, a multicell method was used. The control rod material is placed in the center of eight fuel rods, since this is the configuration that will be used in operation of the NSCR. The cross sections for the control material are homogenized based on the fluxes generated in the fuel. Each shim safety rod is composed of three portions: a 35.56 cm (14 in) long poison region, a 38.10 cm (15 in) long fuel follower under the poison region, and two void regions on the top and bottom with 12.70 cm (5 in) and 17.78 cm (7 in) in length, respectively. The poison in the shim safety rod is borated graphite (B4 C 25%, graphite 75% by weight) and its liameter is 3.3909 cm. The outside clad is 0.0508 cm thick stainless steel. The fuel ! follower is similar to a fuel rod except that the fuel meat is 3.3909 cm in diameter which is smaller than the fuel rod. The void region has the same diameter as the poison region, but is filled with void instead of poison. This geometry is shown in Figure 52. l The transient rod only has two portions: a 38.10 cm (15-in) long poison region and a ; 50.80 cm (20 in) long void region. The poison used in the transient rod is also borated graphite. The poison region has a diameter of 3.03276 cm. The thickness of the stainless steel clad is 0.0712 cm. The geometry of the transient rod is shown if Figure 53. i
l l 90 The regulating rod is composed of a single 38.10 cm (15 in) long poison. The rod size is the same as that of the transient rod. The only difference is that the poison in the ! regulating rod is pure B4C powder (density: 2.5 g/cm3). This geometry is shown in i Figure 54. i i TABLE 13 Atom Densities for the Fuel Rod Model Materials Components Atom Densities (atoms / barn-cm) L Zirconium Rod Zirconium 4.2909E-2 l Uranium-235 8.8382E-4 Uranium-238 - 3.7878E-4 l iP Hydrogen in ZrH 5.4364E-2 fue) f Meat Zirconium 3.3977E-2 Erbium-166 1.0560E-4 : Erbium-167 7.2132E-5 ! Graphite Carbon ' 8.1299E-2 Iron 6.0414E-2 : C Chromium 1.7384E-2 l Fuel Nickel 8.1033E-3 Carbon 3.1716E-4 Iron 5.92%E-2 ! SSW Chromium 1.7385E-2 ; Nickel 7.6995E-3 l Manganese 1.7320E-3 t Aluminum Aluminum 6.0240E-2 ; Hydrogen in H2O 6.6691E-2 Water Oxygen 3.3346E-2 l
.~ ~ , - - - .- . . . - . .-
l 91 Materials Volume Ratio Elements Iron Top Chromium Structure SS304 24.7 % l Nickel of Fuel Manganese Rod Aluminium 17.9 % Aluminum Water 57.4 % Oxygen Rod . n
.h m - ~
C
'l 0.0508em 1r l ! y; l [ M 3.9520cm 5l t 38.10cm - --
3 i l Fuel .
,jjp} : Zirconium Rod J
Rod W 0:45726 min 3r e'~ ~^^7-I Cell Fuel Meat l
' 8 h9C*
n 3, p Ci
.,'- Xg7[
37:vyym, yli 8 3.9520cm F f -- f.4
- WJ f r*. 117en E Ciad I 8.89cm u - . I gb k Water 1r i i 44,iggngg, <
0.050Sem
]
Materials Volume Ratio Elements Bottom Stnicure ' * "
- SS304 32.5 %
Nickel Fuel Rod Manganese Hydrogen in H2O r Water 67.5 % Oxygen Fig. 51. TRIGA FLIP fuel rod geometry used in WIMSD4/m. l
1 i. 1 92
\
1 J ]? 3 1 b Void d L69N5hi. - 3.9520cm g Clad:SS304 .
'(..'
D Water 0.0508cm l J' i n- ~i j 3.9520cm "' ! i l 4 [ E Poison: I 12.70cm , Dorated Graphite l BtC 25%in Wei '
} Graphite 75%
l 3.9520cm N Ciad:SS304 j 35.56cm l 4 0.0508cm W ter l y N 3.9520cm g l 1 i 38.10cm ]? ! & Zirconium Rod v l
$ FuelMeat
, 1r 3 9520cm 1 .. ( A E cia li l7.78cm 0 508cm O w ater ' 1r
-- 3;
- l1 3.9520cm >l f
i Fig. 52. NSCR shim safety control rod geometry used in WIMSD4/m.
?
i 4 l 4 i
. _ - - . . . . - - . . . . . . . . . . . . = . . - - . - . . . . . - - - . - . . - - . . - . - . . - . _ _ - - - . -. .
1 i i i 2 93 s 1 1 .j l 5 I i < I i Poison : l Borateg 9, . i
.. W 25% in EE I't Graphite 75%
3.9520cm g Clad: Aluminum 38.10cm j l j 00712cm i if m Water Ll 7 I- <e tr 3.9520cm -s l .. t
- n 4
l . I I 1 50.80cm '. [
- Void
.li :1:5)638cin f i i
- E Ciad: Aluminum I
I [. ' O water 0.0712cm l
.L l,
! l4 3.9520cm :l l 4 i
- Fig. 53. Transient control rod geometry used in WIMSD4/m. ;
r i 1 I l
94
][ ]I Poison :
Pure B4C Powder Density : 2.5 g/cm3 38.10cm
)b s
9 3.9520cm >l 3.9520cm
][
Clad : Aluminum Water Fig. 54. NSCR regulating control rod geometry used in WIMSD4/m. Other Core Materials The top structure of the fuel rod includes a volume ratio of 24.7% for the stainless steel plug,17.9% for the aluminum handle and 57.4% for the surrounding water. Similarly, in the bottom structure, approximately 32.5% by volume is stainless steel and 67.5% by i volume is water. The size of the parallelepiped graphite reflector elements used on the north and south l sides of the NSCR is 7.62 cm (3 in.) x 7.62 cm (3 in.) x 71.17 cm (28 in.). When placed next to the core, the graphite elements are surrounded by water. Therefore, in the area of the homogenized north and south graphite reflectors 92.94% by volume is graphite and the l rest is water. The graphite was modeled in WIMSD4/m using a multicell model with other graphite cells adjacent to the cell ofinterest on three sides, and fuel adjacent to the cell on l the fourth side. The needed neutronic data for water is obtained by setting up a WIMSD4/m multicell case which is similar to that used to model the control rods. In this case, the central cell is filled by pure water (density: 1.0g/cm ,3 temperature: 300 K) and the other eight cells are fuel. l l
. . . ~ . . . . -. =~.. - - . ..- , ~ ~ _ . _ - - - -
95 APPENDIX H NCTRIGA ; NCTRIGA was used to perform the steady state thermal hydraulic analysis in this thesis. A description of many of the parameters used in this model can be found in Chapter 111. Presented here is a description of the all input parameters and a sample input , deck. A listing of all of the input parameters necessary to run NCTRIGA is shown below. l NCTRIGA Input Parameters Variable Name Definition
- CHANHT Coolant channel height, m.
CHIMMY Unheated section of the reactor above top of fuel j rods, m. CLADTK Thickness of the clad, m. : CPWR Fixed power of whole core, kW. j DELTA Criterion for convergence o f frictional and buoyant forces. DEPTH Distance from pool surface to the bottom of fuel rods, m. l FKVT(I) Fuel conductivity, W/mK, as a function of l temperature, TEMPK(I). 1 FUELHT Length of fuel meat, m. , I l FUELRD Radius of fuel meat, m. FW, FQ, FH Hot channel factors for flow (or enthalpy), heat flux, , and heat transfer coefficient.
96 GAPCON Conductance, W/m2K, between the fuel rod and clad. ICONFG Code for orientation of pins. If = 3, a triangular pitch geometry is assumed; if = 4 a square pitch is assumed. ICTC With a va!ue ofzero the clad conductivity ofIconel 800 is assun.ed and computed internally in the code. Otherwise, the user must input a single, temperature independent value of clad conductivity as UCLADK. IPRT A value of zero denotes an intermediate interation and other output will not be printed. Otherwise, the output will be printed. IQVZ With a value of zero a standard axial power distribution is calculated. Otherwise, the user must input the axial distribution as arrays: QVZ(I) and ZR(I). NK No. of data values in fuel conductivity, FKVT(I), TEMPK(I). With a value of zero the fuel conductivity of Zirconium Hydride is assumed and computed internally in the code. Otherwise, the user must input a temperature dependent array of fuel conductivity as FKVT(I) and TEMPK(I). NN Number of axial nodes in reactor coolant channel, t NTRGEL. Number of TRIGA elements in the whole core. ONEP5 Additive parameter for the code computed fuel conductivity relation that ranges from +1.5 to -1.5. For nominal value use 0.0. POVERD Pitch to diameter ratio (P/D). QVZ(1) Axial power distribution whose averages is one. ZR is the axial average position array. RPEAK Power factor by which this channel differs from an average channel. l l t
?
l 97 l t TEMPK(1) Temprattre array corresponding to FKVT(I). f TPOOL Average temperature of water in pool, C. I UCLADK User input value of clad conductivity, W/mK. ! UHZB, UHZT Unheated lengths at bottom and top of fuel pin, m. j i VGUESS Initial guess for inlet velocity of coolant, m/s. ! I ZR(I) Axial position normalized to one (Z/FUELHT). i Shown below is a sample NCTRIGA input deck. This is included to show specifically , how the input can be constructed to model a channel of the NSCR. i , i I l i J 1 l l I l l i i l i ?- - - - - - - - . ___ _ . - -- , ,., - , , , - . +- , , , .
t i 98 NCTRIGA Input Deck i TAMU NSC TRIGA HEU CORE prepulsing case axial dist from 300 watt-E5 l SW
& INPUT !
?
NN=15 i
, NK=5
, ' NTRGEL=90
, IPRT=0 i
, IQVZ=1 :
, IFK=1 i
, ICTC=1
, ICONFG=4 i
, CPWR=1000.0 !
, TPOOL=37.
, UHZB=0.1778 ;
, UHZT=0.2032 ;
, FUELHT=.381 l'
, FUELRD=.017411
, ONEP5=0 !
, GAPCON=22111 -!
, POVERD=1.10268 !
, CLADTK=.000508 l
, UCLADK=16.8 :
, CHANHT=.5588 l l , CHIMNY=0.000 l
, DELTA =.00001 ;
, . VGUESS=0.2 ;
i
, DEPTH =8.84 5
, RPEAK=1.65778649 ;
L , FW=1.0 l
, FQ=1.0 I
, FH=1.0 FKVT=18.0,18.0,18.0,18.0,18.0, TEMPK= 0.0, 200.,400., 600., 800.,
ZR=0.000, 0.0667, .01333, 0.2000, 0.2667, 0.3333, 0.4000, 0.4667, ; L 0.5333,0.6000,0.6667,0.7333,0.8000,0.8667,0.9333,1.0000, l QVZ=0.4734, 0.5504, 0.661 1, 0.7999, 0.9338, 1.0596, 1.1666, 1.2451, 1.2895,1.2966,1.2652,1.1961,1.0922,0.9609,0.8426,0.8074
&END .
I F --
i l 99 j l APPENDIX C l i PARET Below is a listing of the parameters needed in the construction of a PAIET input deck. Each section contains the values used, along with a description of the methods used to produce these parameters. The headings given below are the same as those used in the i PAIET manual for ease of reference. RESTART CARD The first value is the identifier to indicate a restart problem 1, or initial problem 0. The second value is the frequency with which the restart file should be written in number of time steps. 0 200 Title Line This isjust the name of the job to be printed on the header. The first character must be an asterisk.
- PARET: NSC TRIGA $1.80 pulse at 37 deg C 300 watts II. GeneralInformation This information should be input to the 10xx series cards in the format described above.
l Item Variable Value Description 1 NCHN -2 For the NSCR model, there are two channels. The first is the hot channel with the volume of one actual channel in the core, and the second ! is the average channel having the volume of the rest of the core. The minus signals that SI units will be used. 1 1 l l
3 l 100 j
- 2 NZ 15 In this model, there are 15 axial nodes in the fuel meat. Each is 1 in. ,
3 NR 7 There are seven radial nodes. Five in the fuel, one for the gap, and one for the clad. i l 4 IGEOM 1 Cylindrical geometry . 5 IPROP 1 Reactivity specified problem. If 0 were chosen, the power level could be specified in table 9 instead of the reactivity. 6 IRXSWT 1 This lets the code calculate the quality. 7 IPOP 0 Inlet pressure is specified. 8, KINTS 0 This allows the code to compress the time step, if necessary for kinetics calculations. 9 IDLYGP 6 Six delayed neutron groups. 10 KINPRT 0 No intermediate printout of kinetics parameters. I1 ISUPPR 0 No average temperature printout 12 MAXHCC 20 Limits heat transfer interations to 20 per node. 13 "OWER 0.3-3 Initial Reactor Power in MW. 3 14 PF 0.032656 Total volume of fuelin the core (m ) 15 PRESUR 1.90327+5 Operating pressure in Pa.
-16 . ENTHIN -37 J0 Inlet moderator temperature in -deg C. If a positive value is entered, it represents enthalpy 17 RS 0.01792 Fuel pin radius in m
( l8 RF 0.017411 Radius of fuel meat in m l f
1 I i 101 i 19 RC 0.01742 Radius to inner surface of clad in m ! I 20 PW 0.0 Plate width, doesn't apply here. I i 21 FW 0.0 Used only for plates I i 22 AL 0.381 Active length of fuelin m 23 ALDDIN 0.1778 Length of unfueled inlet, m. 24 ALDEEX 0.2032 Length of unfueled exit, m. l 25 BBEFF 0.0071 Beta effective 26 FL 27.9-6 Prompt neutron generation time, s. 27 GRAV 9.80664 Acceleration due to gravity, m/s2, 28 QW 0.00679 Heat source for moderator. This is the fraction of heat generated in the moderator, multiplied by the ratio of fuel meat volume to i moderator volume. This value gives about l 2% of the energy generation occurring in the moderator. 29 TRANST 15.0 Total time ofinvestigation, s. 30 RXXCON 0.80 A parameter used in void volume generation. 31 RXXEP 1.0 Exponent used in void volume generation. 32 RHOREF 988.67 Moderator reference density kg/m3 . This is the density of the moderator at the initial condition. 33 GAMMA 0 -0.47296 All GAMMAx are coefficients used in fitting the reactivity feedback as a function of temperature. These must be generated from the reactivity changes determined from l DIF3D runs at different temperatures. A detailed description is given in Chapter III. l l
102 34 GAMMA 1 -0.20201-2 35- GAMMA 2 1.15817-5 36 GAMMA 3 0.00 37 GAMMA 4 0.00 38 DOPPN 1.00 Exponent used in reactivity fit. , ! 39 EPS4 0.001 Upperlimit for kinetics time test. l 40. DNBQDP 0.00 This allows the code to calculate the DNB ! heat flux at each node. 41 TAUUNB 0.001 Nucleate boiling bubble collapse time. This value came with the sample deck. 42 TAUUTB -0.001 Transition boiling bubble collapse time. This value came with the sample deck. 43 ALAMNB 0.05 Fraction of clad surface heat flux which is utilized in producing vapor in sub-cooled regien. This value came with the sample deck. 44 ALAMTB 0.05 Same as 43 for transition boiling. This value came with the sample deck. 45 ALAMFB 0.05 Same as 43 for film boiling. 46 HTTCON 1.4 Natural convection heat transfer constant for correlation. l This value came with the sample deck. l . 47 HTTEXP 0.33 Natural Convection heat transfer exponent for l. i original correlation. This value came with the
- sample deck.
1-
'III. AdditionalGeneralInformation This data will be input to the 11lx cards.
103 CARD 1111 la PSUBC 0.049768 This is the total area of flow for all the channels in the core in m2 , 2a FACT 2(1) 1.00 This appears to be some parameter that the code no longer uses. It is recommended to use 1.00 to keep the code from crashing. This value must be repeated for the number of channels used. CARD 1112 lb IONEP 0 The Dittus-Boelter correlation is used for single phase heat conduction 2b- ITWOP 1 The McAdams correlation is used for two phase heat transfer 3b IMODE I The transition modelis used. 4b ICHF 2 The Mirshak DNB correlation is used. 5b IHT 0 The original single phase heat transfer routine is used. 6b QAVE O No value is needed for the average core heat flux. CARD 1113 1c RDRATE 3.81 Rate of control rod movement m/s. This parameter has no effect on the rod movement for the case described by this input deck. The reactivity insertion rate is interpolated from the 9000 cards. 2c TDLAY 0.2 Delay time before control rod movement begins after trip detection, (seconds). 3c POWTP 10000.0 Overpower trip setpoint(MW).
l 1 i 104 )l l l 4c FLOTP 0.00 Low flow trip setpoint. CARD 1114 1d liNCTOP 0.0889 Height above the core considered for natural l convection effects. ' 2d HNCBOT 0.0889 Height below the core considered for natural convection effects. 2000 SERIES CARDS These cards are used to describe the thermal conductivities and heat capacities of the materials used in the fuel. They are fit parameters that must be generated for each type of material. See equations 71 and 72 of the PARET manual for details.' i The following numbers can be used for FLIP fuel. CARD 2001 l l l ALPHA 1, ALPHA 2, ALPHA 3, ALPHA 5 = 0.0 i ALPHA 4 18.00 This is the constant thermal conductivity of the fuel meat (W/mK). General Atomic has j l found that the thernul conductivity of the fuel i is independent of temperature. " i CARD 2002 BETA 1, BETA 4 = 0.00 BETA 2 0.417+4 This is the coefficient of the T term in the fit l equation of the heat capacity. . t i BETA 3 2.04+6 This is the constant in the fit equation of the j i heat capacity. BETA 5 -273.15 This term adjusts the T to units of *C instead , of K which the code generates to match the units of the curve fit. 3
105 CARDS 2003,2004 These are the values used in the fit of the gap conductivity and heat capacity. These are assumed to be constant with respect to temperature. ALPHA 3 0.19900 BETA 3 6.66340+2 CARD 2005,2006 These are the values used in the fit of the clad conductivity and heat capacity. These are assumed to be constant with temperature l ALPHA 3 16.8 BETA 3 3.975+6 3000 CARDS l These cards give the radial description of the fuel pin in the PARET model. For the i TRIGA, seven radial regions are modeled. There are five nodes in the fuel meat, one in the gap, and one in the clad. The first value on each card describes the thickness of the region in meters, the second value denotes the last region with this thickness. The third value denotes the material composition of the region, as described in the 2000 series cards, i.e.1 is fuel,2 is gap,3 is clad. The fourth value l denotes the fraction of the total energy generated in this channel that is generated in i the particular region of the channel. Below are the cards. 3001 is for fuel,3002 for gap,3003 for clad. 3001, 4.35275-3 5 1 0.980 3002, 0.9-5 6 2 0.00 3003, 5.00-4 7 3 0.000 The value of 0.980 on the 3001 card means that 98% of the energy is generated in the fuel. The other 2% is generated in the coolant. This coincides with the value input as item 28, QW, on the 1000 series cards. 4 4000 CARDS
106 These cards are used to describe the axial geometry in the active region of the fuel pin. In this model there are 15 regions of one inch (0.0254 m) each. The card is given below. 4001, 2.54-2 1 2.54-2 14 2.54-2 15 5000 CARDS These cards give information about each channel. For the first channel, the cards are numbered from 5100, for the second channel, the cards are numbered from 5200. The data is input as follows. I i ITEM VARIABLE CH1 VALUE CH2 VALUE DESCRIPTION 2a IFLOW l 1 The flow for all channels will be specified on the 10000 series cards as a function of time. 3a DELP 0 0 Not needed for this model 4a RN 0.02794 0.02794 This is the distance in meters from l the center of the pin to the node m ; the center of the water channel. I For pin geometry, this is taken to ! be pitch /2*sgrt(2) Sa BM 1.111-2 9.999-1 This is the reactivity weighting ; factor for the channel. These ' values make 1/90 of the core from ; the hot channel, and 89/90 of the , core from the average channel.
]
6a ALOSCN 0.5 0.5 Unrecoverable loss coefficient for the inlet. 7a ALOSCX 0.55 0.55 Unrecoverable loss coefficient at the exit. 8a SIGIN 1.00 1.00 Ratio of channel area to inlet plenum area.
_ _ . _ . ._ . ~ _ _ . _ _ _ . _ _ _ _ . _ _ _ _ . _ _ . _ _ _ _ l I- ! 107 l l l l l 9a SIGEX 1.00 1.00 Ratio of channel area to outlet plenum area. l 10a DVOID 0.00 0.00 Void / density coefficient that is not needed here. j 11a DTMP 0.00 0.00 Coolant temperature coefficient ! that is not used in this model. 1 CARD 5101,5201 l l These cards represent the plenea. These same values are used for both channels ' and are as follows. , Ib ALPPIN 0.0889 This is the length of the inlet plenum. 2b ALPPEX 0.0889 This is the length of the outlet plenum. ! I s 3b DEEIN 4.38511-2 This is the equivalent diameter of the inlet plenum. ! 4b DEEEX 4.38511-2 This is the equivalent diameter of the outlet f plenum. : I The remaining 5000 series cards contain information describing the axial power description, followed by weighting factors which have all been set to 1.00 since , they are not used in this model. The axial power description is the ratio of the local ! to average power. This is not normalized to one, except for the average channel. i The hot channel should be normalized to the radial peaking factor. Fifteen values j must be input. The first and the last are the peaking factors for the top and bottom , of the fuel, not the nodes output from DIF3D. These values must be extrapolated - from data generated by DIF3D. The position of the other nodes in PARET corresponds exactly to the nodes output from DIF3D, so no interpolation is l necessary. The value input on these cards represent the state of the reactor prior to ; the transient.. ; L The entire set of 5x00 cards must be repeated for each channel modeled. A complete set of 5100 cards and 5200 cards are shown below. 5100, .1 0 0.02794 1.1111-2 0.5 0.55
- 5100, 1.00 1.00 0.00 0.00 ;
5101, 0.0889 0.0889- 4.38511-2 A38511-2 i I
( l 108 5102, 0.6563 1.00 1.00 1.00 5103, 0.8587 1.00 1.00 1.00 5104, 1.0636 1.00 1.00 1.00 5105, 1.2791 1.00 1.00 1.00 5106, 1.5037 1.00 1.00 1.00 5107, 1.6854 1.00 1.00 1.00 5108, 1.8191 1.00 1.00 1.00 5109, 1.9022 1.00 1.00 1.00 5110, 1.9326 1.00 1.00 1.00 5111, 1.9091 1.00 1.00 1.00 5112, 1.8314 1.00 1.00 1.00 5113, 1.7007 1.00 1.00 1.00 5114, 1.5208 1.00 1.00 1.00 5115, 1.3093 1.00 1.00 1.00 5116, 1.1314 1.00 1.00 1.00 5200, 1 0 0.02794 9.9999-1 0.5 0.55 5200, 1.00 1.00 0.00 0.00 5201, 0.0889 0.0889 4.38511-2 4.38511-2 5202, 0.5087 1.00 1.00 1.00 5203, 0.6323 1.00 1.00 1.00 5204, 0.7771 1.00 1.00 1.00 5205, 0.9084 1.00 1.00 1.00 5206, 1.0343 1.00 1.00 1.00 5207, 1.1393 1.00 1.00 1.00 l 5208, 1.2172 1.00 1.00 1.00 5209, 1.2641 1.00 1.00 1.00 5210, 1.2776 1.00 1.00 1.00 5211, 1.2567 1.00 1.00 1.00 5212, 1.2009 1.00 1.00 1.00 5213, 1.1112 1.00 1.00 1.00 5214, 0.9902 1.00 1.00 1.00 5215, 0.8505 1.00 1.00 1.00 5216, 0.7438 1.00 1.00 1.00 6000 CARDS i These cards provide the information concerning the delayed neutron emission. The parameter IDLYGP, the ninth item on the 1000 series cards, specifies the number of delayed neutron groups that are used. This data is in ordered pairs; the first 4 number is the yield fraction for the group and the second number is the decay constant. 4
109 This model contains six groups. The data is given below. 6001, 0.3824-1 0.12722-1 0.21194 0.31737-1 0.1878 0.1161 6002, 0.40684 0.31137 0.12879 1.4001 0.2639-1 3.8706 9000 CARDS These cards represent the reactivity as a function of time in this model. They can also represent the power level as a function of time, depending on the type of problem desired. The 9000 card contains one integer value signaling the number of i sets of data contained on the subsequent cards. The minimum number of data sets is two. The remaining cards begin their numbering at 9001. The first value is the reactivity desired and the second value is the time at which this reactivity will be reached. PAIET performs linear interpolation for times between the input values. The first value input must be for time 0.00 and the last time must be greater than the total time of the transient. Sample cards for a $1.80 pulse are shown below. 9000, 6 9001, 0.000 0.000 0.000 0.15 -1.80 0.20 9002, 1.800 2.000 0.00 2.20 0.00 1000.0 Note that the reactivity is inserted between 0.15 and 0.20 seconds and is removed between 2.000 and 2.20 seconds. This represents the transient rod ejection under pressure, which causes faster rod motion than the insertion under the influence of gravity. 10000 CARDS These cards represent the inlet mass flow rate as a function of time. These cards work much like the 9000 cards. Here the 10000 card consists of one integer value denoting the number of ordered pairs on the subsequent cards. The cards numbered from 10001 consist of ordered pairs of the flow rate,in units of kg/m 2s, and the time in seconds at which the flow rate takes effect. Again, the code interpolates linearly between the values for intermediate times. Sample cards are shown below. i l
I i 110 l l 10000, 6 10001, 10.49 0.00 300.0 0.30 10002,' 300.0 0.50 100.0 5.00 10003, 10.49 45.0 10.49 1000.0 11000 CARDS These cards represent the thermal expansion of the clad as a function of temperature. This is another series of cards with the first card telling how many sets of data are contained on the other cards, and the next cards containing data. The data is entered in pairs, the first parameter is the percent of expansion of the clad. The second is the temperature in *C at which this percent expansion exist. Sample cards are shown below. I1000, 2 11001, 0.0 10.00 0.5 2000.0 12000 CARDS These cards are not used in this model. The pressure drop as a function of time is calculated by the code. All values are set to z.ero. 14000 CARDS These cards represent the calculational time step as a function of transient time. This is important for producing accurate results with optimum efficiency. A small time step is desired when the power or temperature is changing rapidly, and a large time step is desired when the parameters are changing slowly. Sample cards are shown below. 14000, 6 14001, 0.001 0.0 0.0001 0.20 0.00005 0.25 14002, 0.001 1.0 0.001 2.20 0.005 3.0
- 16000 CARDS These cards determine the frequency which the major and minor edits are output. The values input here are determined by how much resolution is desired in the printed results. This does not affect the calculations; the code will continue to l
l
l 111 calculate in the steps given above, but it will not print the results at every time step. The output files from PARET are typically quite large. The size of the file produced using the following cards for the 15 second transient is about 15 Mb. The l cards work in sets of three numbers. The first is the time increment for major edits, the second is the frequency on intermediate edits, and the third is the time at which these go into effect. i 1 16000, 5 , 16001, 0.10 50 0.0 0.0010 10 0.10 16002, 0.05 10 1.00 0.005 20 2.50 16003, 0.5 10 3.0 17000 CARDS These cards represent relative pump mass velocity as a function of time. For I the NSCR case, there is no forced convection. All values are set to 1.00. ) I 17000, 2 l 17001, 1.0 0.0 1.00000 455.0 l 18000 CARDS These card.i represent the worth of the control elements as a function of position. These cards have no effect on the results when the reactivity is input as a function of time using the 9000 cards. 18000, 2 l 18001, 0.0 0.0 -4.00 0.381 l l t__________ _ _ _ _ _
I i i i 112 l l . l APPENDIX D ; l PROCESSING CODES ! i Two data processing codes were written to facilitate repetitive processes in the i analyses. The first is Power 2, which reads output edits from DIF3D and generates power I peaking factors averaged over a fuel rod. The second code is AXPOW, which reads the l output from Power 2 and writes input cards for PARET and NCTRIGA. Both Power 2 and j j AXPOW were written in FORTRAN and compiled and executed on a Sun SparcStation. ( Power 2 POWER 2 was written primarily by Chien-Hsaing Chen. It reads the three dimensional power distribution edit from DIF3D and writes averaged power distributions axially and radially for each fuel cell. The source listing for this code is shown below. . CHARACTER XCOOR(5)*1,MAXDIR*2,maxcord*1 REAL*8 V(120,100,100),
+ AP NW(10,10), AP N E(10,10),APSW(10,10), AP S E(10,10),
+ PDNW,PDNE,PDSW,PDSE,PD,
+ NW(10,10,15),NE(10,10,15),SW(10,10,15),SE(10,10,15),
+ RPNW(10,10,15),RPNE(10,10,15),
+ RPSW(10,10,15),RPSE(10,10,15) ,
INTEGER I,J,K,lP ; O P E N (10, FI LE =' power. in p',STATU S='O LD') I OPEN(20, FILE =' power.dat', STATUS ='NEW) ; C ! C READ POWER C ; DO K=1,51 : DO IP=1,4 READ (10,100) DO J=49,1,-1 READ (10,200)(V(1,J,K),1=(IP-1)*9+1,lP*9) ENDDO ENDDO ' READ (10,100) DO J=49,1 -1 READ (10,200)(V(1,J,K),1=37,40) i ENDDO l
i 113 ; ENDDO ! C DO K=1,15 l NW(1,7, K) = .25*V(7,10, K+ 18) +.25 *V(8,10, K+ 18) +
+ .25*V(7,1 1, K+ 18)+ .25*V(8,1 1, K+ 18) i N E(1,7, K) =.25*V(9,10, K+ 18) +.25*V( 10,10, K+ 18) +
+ .25*V(9,1 1, K + 18)+.25*V(10,1 1, K+ 18) ,
SW(1,7,K)=.25N(7,12,K+18)+.25*V(8,12,K+18)+ i'
+ .25*V(7,13,K+18)+.25*V(8,13,K+18)
S E(1,7, K) =.25*V(9,12, K+ 18)+.25 *V( 10,12, K+ 18) +
+ .25*V(9,13, K+ 18) +.25*V(10,13, K+ 18)
NW(2,7, K) =.125*V( 1 1,10, K+ 18)+.125*V( 12,10, K+ 18)+
+ .125*V(13,10, K+ 18) +.125N(14,10,K+ 18) + !
i
+ .125*V( 1 1,1 1, K+ 18) +.125*V(12,1 1,K+ 18) +
+ .125*V(13,11,K+18)+.125*V(14,11,K+18) l N E(2,7,K)=.25*V(15,10,K+ 18)+.25*V(16,10,K+ 18)+
+ .25*V(15,11,K+18)+.25N(16,11,K+18) !
SW(2,7, K)= .125*V(1 1,12, K+ 18)+.125*V(12,12,K+ 18)+
+ .125*V(13,12, K+ 18) + .125*V(14,12, K+ 18)+
+ .125*V(11,13,K+ 18)+.125N(12,13, K+ 18) +
+ .125*V(13,13, K+ 18) +.125N(14,13, K+ 18)
S E(2,7, K) =.25*V(15,12,K+ 18)+.25*V(16,12, K+ 18) +
+ .25*V(15,13, K+ 18) +.25 *V(16,13,K+ 18)
NW(4,7, K)=.125 *V(23,10, K+ 18) +.125*V(24,10, K+ 18) +
+ .125*V(25,10, K+ 18)+.125*V(26,10,K+ 18) + ;
+ .125 *V(23,1 1,K+ 18) + .125*V(24,1 1,K+ 18) +
+ .125 *V(25,11, K+ 18) +.125*V(26,11, K+ 18)
N E(4,7, K)= .25 *V(27.10, K+ 18) + .25*V(28,10,K+ 18) +
+ .25 *V(27,1 1, K+ 18) +.25*V(28,1 1, K+ 18) ;
SW(4,7, K)=.125*V(23,12, K+ 18)+.125*V(24,12, K+ 18) +
+ .125*V(25,12,K+ 18)+.125*V(26,12,K+ 18)+
+ .125*V(23,13, K+ 18) +.125*V(24,13,K+ 18) +
+ .125*V(25,13, K+ 18) +.125*V(26,13, K+ 18)
S E(4,7,K)=.25*V(27,12,K+ 18)+.25*V(28,12,K+ 18)+
+ .25 *V(27,13, K+ 18)+.25*V(28,13, K+ 18)
NW(5,7,K)=.125*V(29,10,K+18)+.125*V(30,10,K+18)+
+ .125*V(31,10, K+ 18)+.125*V(32,10, K+ 18) +
+ .125*V(29,1 1, K+ 18) + .125*V(30,1 1, K+ 18) +
+ .125*V(31,1 1, K+ 18)+.125*V(32,11, K+ 18)
N E(5,7, K)=.25 *V(33,10, K+ 18) +.25*V(34,10,K+ 18) + l t + .25*V(33,11, K+ 18) +.25*V(34,11,K+ 18) SW(5,7,K)=.125*V(29,12,K+18)+.125*V(30,12,K+18)+ l + .125*V(31,12,K+18)+.125*V(32,12,K+18)+ l + .125*V(29,13,K+ 18)+ .125*V(30,13,K+ 18)+
- j. + .125*V(31,13,K+18)+.125*V(32,13,K+18)
! S E(5,7,K)= .25N(33,12, K+ 18) +.25*V(34,12,K+ 18)+
+ .25 *V(33,13, K+ 18)+.25*V(34,13, K+ 18) l N E(1,6, K) =.25*V(9,14, K+ 18)+.25*V(10,14,K+ 18)+
l + .25*V(9,15,K+ 18) +.25*V(10,15,K+ 18) S E(1,6,K)=.125*V(9,16, K+ 18) + .125*V(10,16,K+ 18) + l + .125*V(9,17,K+18)+.125*V(10,17,K+18)+ ! + .125N(9,18, K+ 18) +.125*V(10,18,K+ 18) +
+ .125*V(9,19,K+ 18)+,125*V(10,19,K+ 18)
- NW(2,6, K)=.125*V(1 1,14, K+ 18) +.125*V(12,14,K+ 18) +
i + .125*V(13,14,K+ 18)+.125*V(14,14, K+ 18) + j + .125*V(1 1,15, K+ 18)+.125*V(12,15,K+ 18) +
+ .125*V(13,15,K+ 18) +.125*V(14,15, K+ 18) i r
114 i N E(2,6,K)=.25*V(15,14, K+ 18)+.25*V(16,14,K+ 18)+
+ .25*V(15,15, K + 18) +.25 N(16,15, K+ 18)
SW(2,6, K)= . 0625 *V( 1 1,16, K+ 18) + . 0625*V( 12,16, K+ 18) +
+ .0625*V(13,16,K+ 18) +.0625*V(14,16,K+ 18)+
+ .0625*V(11,17,K+ 18)+.0625*V(12,17,K+ 18)+
+ .0625 *V(13,17, K+ 18) + .0625N(14,17, K+ 18) +
+ .0625 N(11,18, K+ 18) +.0625 *V(12,18, K+ 18)+
+ .0625N(13,18, K+ 18) +.0625 *V(14,18, K+ 18) +
+ .0625*V(11,19, K+ 18) +.0625 *V(12,19, K+ 18)+
+ .0625*V(13,19, K+ 18)+.0625*V(14,19, K+ 18)
S E(2,6, K)=.125 *V( 15,16, K+ 18)+ .125*V(16,16, K+ 18) +
+ .125 *V( 15,17,K+ 18)+ .125 *V( 16,17, K+ 18)+
+ .125*V(15,18, K+ 18) + .125 *V(16,18, K+ 18)+
+ .125 *V(15,19, K+ 18) +.125*V(16,19, K+ 18)
NW(3,6, K) =.125*V( 17,14, K+ 18) +.125*V( 18,14, K+ 18) +
+ .125*V(19,14, K+ 18)+ .125*V(20,14, K+ 18)+
+ .125*V(17,15,K+ 18)+.125*V(18,15,K+ 18)+
+ .125*V( 19,15, K+ 18) +.125*V(20,15, K+ 18)
N E(3,6,K)=.25*V(21,14, K+ 18)+.25*V(22,14,K+ 18)+
+ .25*V(21,15, K+ 18) +.25*V(22,15, K+ 18)
SW(3,6,K)=.0625*V(17,16, K+ 18)+.0625*V(18,16,K+18)+
+ .0625 *V(19,16, K+ 18) +.0625 *V(20,16, K+ 18)+
+ .0625 *V(17,17, K+ 18)+.0625 *V(18,17, K+ 18)+
+ .0625 N(19,17, K+ 18) +.0625 N(20,17, K+ 18) +
+ .0625*V(17,18, K+ 18) +.0625*V(18,18, K+ 18) +
+ .0625*V(19,18, K+ 18) +.0625 *V(20,18, K+ 18) +
+ . 0625*V(17,19, K+ 18) + .0625 *V( 18,19, K+ 18) +
+ .0625*V(19,19, K+ 18) +.0625 *V(20,19, K+ 18)
S E(3,6, K)=.125*V(21,16, K+ 18)+.125*V(22,16, K+ 18)+
+ .125*V(21,17, K+ 18)+.125 *V(22,17, K+ 18) +
+ .125*V(21,18,K+ 18)+.125*V(22,18,K+ 18)+
+ .125*V(21,19, K+ 18) + .125*V(22,19, K+ 18)
NW(4,6, K) =.125*V(23,14, K+ 18) + .125*V(24,14, K+ 18) +
+ .125*V(25,14, K+ 18)+.125*V(26,14, K+ 18) +
+ .125*V(23,15, K+ 18) + .125*V(24,15, K+ 18) +
+ .125 *V(25,15, K+ 18) +.125*V(26,15, K+ 18)
N E(4,6,K)=.25*V(27,14,K+ 18)+.25*V(28,14,K+ 18)+
+ .25*V(27,15, K+ 18) +.25N(28,15, K+ 18)
SW(4,6,K)=.0625*V(23,16,K+18)+.0625*V(24,16,K+18)+
+ .0625*V(25,16, K+ 18)+.0625*V(26,16, K+ 18) +
+ . 0625*V(23,17, K+ 18) +.0625*V(24,17, K+ 18) +
+ .0625N(25,17,K+ 18) +.0625 *V(26,17, K+ 18)+
+ .0625*V(23,18,K+ 18)+.0625*V(24,18, K+ 18)+
+ .0625N(25,18, K+ 18) +.0625*V(26,18, K+ 18) +
+ .0625*V(23,19,K+18)+.0625N(24,19,K+18)+
+ .0625*V(25,19, K+ 18) +.0625 *V(26,19, K+ 18) !
S E(4,6, K)= .125*V(27,16, K+ 18)+.125 *V(28,16, K+ 18) + '
+ .125*V(27,17, K+ 18)+.125 *V(28,17, K+ 18)+
+ .125*V(27,18, K+ 18) + .125*V(28,18, K+ 18)+
+ .125 *V(27,19, K+ 18) +.125*V(28,19, K+ 18)
NW(5,6,K)=.125*V(29,14,K+18)+.125*V(30,14,K+18)+ l + .125*V(31,14, K+ 18) +.125*V(32,14, K+ 18)+ 1
+ .125*V(29,15,K+ 18)+.125*V(30,15,K+ 18) +
i '
+ .125*V(31,15,K+ 18) +.125*V(32,15, K+ 18) 1 N E(5,6,K)=.25*V(33,14,K+ 18)+.25N(34.14,K+ 18)+
+ .25*V(33,15,K+ 18) +.25*V(34,15,K4 18)
SW(5,6, K)=.0625*V(29,16, K+ 18)+.0625N(30,16,K+ 18)+ l I i i
115
+ .0625*V(31,16,K+ 18)+.0625 *V(32,16,K+ 18)+
+ _
.0625*V(29,17,K+18)+ 0625*V(30,17,K+18)+
+ .0625*V(31,17,K+ 18)+.0625*V(32,17,K+ 18)+
+ .0625*V(29,18, K+ 18) + . 0625 *V(30,18, K+ 18) +
+ .0625*V(31,18, K+ 18)+.0625*V(32,18,K+ 18)+ i
+ .0625*V(29,19, K+ 18) + . 0625*V(30,19, K+ 18)+
+ .0625*V(31,19, K+ 18)+.0625 *V(32,19,K+ 18)
S E(5,6, K) =.125*V(33,16, K+ 18) +.125 *V(34,16, K+ 18)+
+ .125*V(33,17, K+ 18) + .125*V(34,17, K+ 18) +
+ .125*V(33,18,K+18)+.125N(34,18,K+18)+
+ .125*V(33,19,K+ 18)+.125*V(34,19, K+ 18)
NW(1,5, K)=.25*V(7,20, K+ 18)+.25 *V(8,20, K+ 18)+ ! + .25*V(7,21, K+ 18) +.25*V(8,21, K+ 18) l N E(1,5, K)=. 25*V(9,20,K+ 18)+.25*V( 10,20, K+ 18)+
+ .25*V(9,21, K+ 18) +.25*V(10,21, K+ 18) l SW(1,5, K)=.125*V(7,22, K+ 18)+.125*V(8,22, K+ 18)+
+ .125*V(7,23, K+ 18)+.125*V(8,23, K+ 18) + .
+ .125*V(7,24,K+18)+.125*V(8,24,K+18)+ l
+ .125*V(7,25, K+ 18) +.125 *V(8,25, K+ 18)
S E(1,5, K) =.125*V(9,22, K+ 18) +.125 *V( 10,22, K+ 18) +
+ .125 *V(9,23, K+ 18) + .125*V(10,23,K+ 18) +
l
+ .125*V(9,24,K+ 18)+.125*V(10,24,K+ 18)+ !
l 1 + .125 *V(9,25, K+ 18) +.125 *V(10,25, K+ 18) i NW(2,5, K) =.125*V( 1 1,20,K+ 18) + .125*V(12,20, K+ 18) + l
- + .125*V( 13,20, K+ 18)+.125*V( 14,20, K+ 18) +
l + .125*V( 1 1,21,K+ 18)+.125 *V( 12,21,K+ 18) +
+ .125 N(13,21, K+ 18) +.125*V(14,21, K+ 18) l N E(2,5,K)=.25*V(15,20,K+ 18)+.25*V(16,20, K+ 18)+
+ .25*V( 15,21, K+ 18) +.25*V( 16,21, K+ 18) l SW(2,5, K)=. 0625*V( 1 1,22, K+ 18) +. 0625*V(12,22,K+ 18) +
i
+ .0625*V(13,22, K+ 18) +.0625 *V(14,22,K+ 18) +
! + .0625*V(11,23,K+18)+.0625*V(12,23,K+18)+ -
+ .0625 *V(13,23, K+ 18) +.0625 *V(14,23,K+ 18)+
+ .0625*V(11,24,K+18)+.0625 V(12,24,K+18)+
+ .0625*V(13,24, K+ 18) +.0625 *V(14,24,K+ 18)+
+ .0625*V(11,25, K+ 18)+.0625 *V(12,25, K + 18)+
+ .0625*V(13,25,K+18)+.0625*V(14,25,K+18) l S E(2,5, K)=.125*V(15,22, K+ 18)+.125*V(16,22, K+ 18)+
+ .125*V(15,23,K+18)+.125*V(16,23,K+18)+
+ .125*V(15,24,K+18)+.125*V(16,24,K+18)+
l
+ .125*V(15,25, K+ 18)+.125*V(16,25, K+ 18)
NW(3,5, K)=.125*V(17,20, K+ 18)+.125 *V(18,20,K+ 18)+ l
+ .125*V(19,20,K+ 18)+.125*V(20,20,K+ 18)+
+ .125*V(17,21,K+18)+.125*V(18,21,K+18)+
+ .125N(19,21, K+ 18)+.125N(20,21, K+ 18)
N E(3,5,K)=.25*V(21,20, K+ 18) +.25N(22,20, K+ 18) +
+ .25N(21,21,K+18)+.25*V(22,21,K+18)
SW(3,5,K)=.0625*V(17,22,K+18)+.0625*V(18,22,K+18)+
+ .0625*V(19,22,K+18)+.0625*V(20,22,K+18)+
+ .0625N(17,23, K+ 18) +.0625N(18,23, K+ 18) +
+ .0625*V(19,23,K+18)+.0625*V(20,23,K+18)+
. + .0625*V( 17,24, K+ 18) +.0625 *V(18,24, K+ 18) + ! + .0625*V(19,24, K+ 18) +.0625*V(20,24,K+ 18)+ ' ' + .0625*V(17,25,K+ 18) +.0625*V(18,25,K+ 18)+
+ .0625*V(19,25, K+ 18)+.0625 N(20,25,K+ 18)
SE(3,5,K)=.125*V(21,22,K+18)+ 125*V(22,22,K+18)+
+ .125*V(21,23, K+ 18)+.125N(22,23,K+ 18) +
l
-- . ~ . . .. . . . . _. - .. - .-
116 r
+ .125 *V(21,24, K + 18) + .125*V(22,24, K+ 18) +
+ .125 *V(21,25,K+ 18)+.125*V(22,25, K+ 18)
N W(4,5, K)=.125 *V(23,20, K+ 18) +.125*V(24,20, K+ 18) +
+ .125*V(25,20,K+ 18)+.125*V(26,20,K+ 18)+ '
+ .125*V(23,21, K+ 18)+.125 *V(24,21, K+ 18) +
+ .125 *V(25,21, K+ 18) +.125 *V(26,21, K+ 18)
N E(4,5, K) =.25*V(27,20, K+ 18) +.25*V(28,20, K+ 18) + 4
+ . 25*V(27,21,K+ 18)+.25*V(28,21,K+18) i SW(4,5, K) =.0625*V(23,22, K+ 18)+.0625*V(24,22, K+ 18) +
+ .0625*V(25,22, K+ 18) +.0625*V(26,22, K+ 18) +
+ .0625*V(23,23, K+ 18) +.0625*V(24,23, K+ 18) +
+ .0625*V(25,23, K+ 18)+.0625*V(26,23, K+ 18) +
+ .0625*V(23,24,K+ 18) + .0625*V(24,24, K+ 18) + !
+ .0625*V(25,24,K+18)+ 0625*V(26,24,K+18)+
+ .0625 *V(23,25, K+ 18)+.0625N(24,25, K+ 18) +
+ .0625 *V(25,25, K+ 18) +.0625'V(26,25, K+ 18) ,
S E(4,5,K)=.125*V(27,22,K+ 18)+.125*V(28,22,K+ 18)+
+ .125 *V(27,23, K+ 18)+.125*V(28,23, K+ 18)+
+ .125*V(27,24, K+ 18)+ .125*V(28,24, K+ 18) +
+ .125 *V(27,25,K+ 18) +.125 *V(28,25, K+ 18)
NW(5,5, K) =.125*V(29,20, K+ 18)+.125 *V(30,20, K+ 18) + ;
+ .125*V(31,20, K+ 18) +.125*V(32,20, K+ 18)+
+ .125*V(29,21, K+ 18)+.125*V(30,21, K+ 18)+ ,
+ .125*V( 31,21,K+ 18) +.125*V(32,21, K+ 18)
N E(5,5,K)=.25*V(33,20,K+ 18)+.25*V(34,20,K+ 18)+
+ .25 *V(33,21, K+ 18)+. 25*V( 34,21, K+ 18)
SW(5,5, K)= .0625*V(29,22,K+ 18) +.0625*V(30,22, K+ 18) + ,
+ .0625*V(31,22, K+ 18)+ .0625 *V(32,22, K+ 18)+
+ . 0625 *V(29,23, K+ 18) +. 0625 *V(30,23, K+ 18)+
+ .0625*V(31,23,K+ 18) + . 0625*V(32,23,K+ 18) + :
+ .0625*V(29,24, K+ 18) + .0625 *V(30,24, K+ 18) +
+ .0625 *V(31,24, K+ 18)+ .0625 *V(32,24, K+ 18) +
+ .0625*V(29,25, K+ 18)+.0625 *V(30,25, K+ 18) + 1
+ .0625*V(31,25, K+ 18) +.0625*V(32,25, K+ 18)
S E(5,5,K) = .125*V(33,22, K+ 18) + .125 *V(34,22,K+ 18) + ,
+ .125*V(33,23,K+ 18)+.125*V(34,23,K+ 18) +
+ .125*V(33,24,K+18)+.125*V(34,24,K+18)+
+ .125N(33,25,K+ 18) +.125*V(34,25, K+ 18)
N E( 1,4, K)=. 25*V(9,26, K+ 18) +.25*V( 10,26, K+ 18) + l
+ .25N(9,27,K+ 18) +.25*V(10,27,K+ 18) i i SE(1,4,K)=.125*V(9,28,K+18)+.125*V(10,28,K+18)+
+ .125*V(9,29,K+ 18) +.125N( 10,29, K+ 18) +
+ .125*V(9,30, K+ 18)+.125*V( 10,30, K+ 18) +
l + .125N(9,31,K+18)+.125*V(10,31,K+18) NW(2,4, K)=.125*V(1 1,26, K+ 18) +.125*V(12,26, K+ 18) +
+ .125*V(13,26, K+ 18)+.125 *V(14,26,K+ 18) +
+ .125*V( 1 1,27,K+ 18) +.125*V(12,27, K+ 18) +
+ .125*V(13,27, K+ 18) +.125*V(14,27, K+ 18)
N E(2,4,K)=.25N(15,26,K+ 18)+.25N(16,26,K+ 18)+
+ .25 *V(15,27, K+ 18) +.25*V(16,27,K+ 18) j' SW(2,4, K) =.0625*V(1 1,28, K+ 18) +.0625*V(12,28, K+ 18) +
+ .0625N(13,28,K+ 18)+.0625'V(14,28, K+ 18)+
+ .0625N(11,29,K+18)+.0625*V(12,29,K+18)+
+ .0625N(13,29,K+ 18)+.0625N(14,29,K+ 18)+
+ 06255'(11,30,K+ 18) +.0625*V(12,30, K+ 18)+
+ .0625*V(13,30,K+18)+.0625N(14,30,K+18)+
+ .0625N(11,31,K+ 18)+.0625 *V(12,31, K+ 18)+
l
_ _ ~. .~ _ _ _. - . _. . _ . 117
+ .0625*V(13,31,K+ 18)+.0625*V(14,31, K+ 18)
S E(2,4,K)=,125*V(15,28, K+ 18)+.125*V(16,28.K+ 18)+ .
+ .125*V(15,29,K+ 18)+.125*V(16,29,K+ 18) + + .125 *V( 15,30, K+ 18) + .125*V( 16,30, K+ 18) + + .125*V(15,31,K+ 18)+.125*V(16,31, K+ 18)
N W(3,4, K) =.125*V(17,26,K+ 18) +.125*V( 18,26, K+ 18)+
+ ,125*V(19,26,K+ 18)+.125*V(20,26,K+ 18)+
l + .125*V( 17,27, K+ 18)+.125 N(18,27, K+ 18)+
+ .125*V(19,27, K+ 18) +.125*V(20,27, K+ 18) l N E(3,4,K)=.25*V(21,26,K+ 18)+.25*V(22,26, K+ 18)+ + .25*V(21,27, K+ 18)+.25*V(22,27,K+ 18)
- SW(3,4,K)=.0625*V(17,28, K+ 18) +.0625*V(18,28, K+ 18) +
l + .0625*V(19,28,K+ 18)+.0625*V(20,28,K+ 18)+
+ .0625N(17,29,K+ 18) +.0625 *V(18,29, K+ 18) +
+ .0625*V(19,29,K+ 18)+.0625*V(20,29,K+ 18)+
+ .0625*V(17,30, K+ 18) +.0625*V(18,30, K+ 18) +
+ .0625N(19,30, K+ 18) +.0625*V(20,30, K+ 18)+
+ .0625*V(17,31, K+ 18) +.0625*V(18.31,K+ 18) + j
+ .0625*V(19,31, K+ 18)+.0625*V{20,31, K+ 18) i SE(3,4,K)=.125*V(21,28,K+18)+ 125*V(22,28,K+18)+
+ .125*V(21,29, K+ 18)+.125*V(22,29, K+ 18)+
+ .125 *V(21,30, K+ 18)+ .125 *V(22,30, K+ 18) + i
+ .125"/(21,31, K+ 19) +.125 *V(22,31, K+ 18) l NW(4,4,K) s.125*V(23,20,K+18)+.125N(24,26,K+18)+ l
+ .125* '(25,26,K+1f.)+.125*V(26,26,K+18)+ 1
+ .125 *VG? ? 7,518)+.125 *V(24,27, K+ 18) +
+ .125 *V(25,27, K+ 18)+.125*V(26,27, K+ 18)
N E(4,4, K)=.25 *V(27,26, K+ 18) + .25*V(28,26,K+ 18) +
+ .25 *V(27,27, K+ 18) +.25 *V(28,27, K+ 18)
SW(4,4, K) =.0625*V(23,28, K+ 18) +.0625 *V(24,28, K+ 18) +
+ .0625*V(25,28,K+ 18)+.0625*V(26,28, K+ 18)+ + .0625*V(23,29,K+ 18)+.0625*V(24,29,K+ 18)+ + .0625*V(25,29, K+ 18) +.0625*V(26,29, K+ 18) + + .0625*V(23,30,K+ 18)+.0625*V(24,30,K+ 18)+ + .0625 *V(25,30,K+ 18) +.0625*V(26,30, K+ 18) + + .0625*V(23,31, K+ 18)+.0625*V(24,31, K+ 18) + + .0625*V(25,31,K+18)+.0625*V(26,31,K+18) -
S E(4,4,K)=.125*V(27,28,K+ 18)+.125*V(28,28,K+ 18)+
+ .125*V(27,29,K+ 18)+.125*V(28,29, K+ 18)+ + .125N(27,30, K+ 18)+.125*V(28,30, K+ 18) + + .125*V(27,31, K+ 18)+.125*V(28,31, K+ 18) !
l NW(5,4,K) =.125*V(29,26,K+ 18) +.125N(30,26,Kv i 8) +
+ .125*V(31,26, K+ 18)+ .125 *V(32,26, K+ 18) + l + .125*V(29,27, K+ 18)+,125N(30,27, K+ 18)+ + .125 *V(31,27,K+ 18)+.125*V(32,27, K+ 18)
N E(5,4,K)=.25*V(33,26, K+ 18) +.25*V(34,26, K+ 18) +
+ .25*V(33,27,K+ 18) +.25N(34,27, K+ 18) i SW(5,4,K)=.0625*V(29,28,K+18)+.0625*V(30,28,K+18)+
+ 0625N(31,28,K+18)+ 0625*V(32,28,K+18)+
+ .0625*V(29,29,K+18)+.0625*V(30,29,K+18)+
+ .0625*V(31,29,K+ 18)+.0625*V(32,29, K+ 18) + ,
j + .0625*V(29,30, K+ 18) +.0625 *V(30,30,K+ 18) + '
+ .0625*V(31,30,K+ 18)+.0625*V(32,30, K+ 18)+
- + .0625*V(29,31, K+ 18) +.0625*V(30,31, K+ 18)+
+ .0625 *V(31,31,K+ 18)+.0625*V(32,31, K+ 18)
S E(5,4,K) =.125*V(33,28,K+ 18) +.125*V(34,28, K+ 18)+
+ .125*V(33,29,K+18)+.125*V(34,29,K+18)+
l l 4 i
-~ _ .- .- - -. . -. . - .. -
118
+ .125*V(33,30, K+ 18)+.125 *V(34,30, K+ 18) +
+ .125 *V(33,31, K+ 18) +.125 *V(34,31, K+ 18)
NW( 1,3, K)=.125 *V(7,32, K+ 18)+.125*V(8,32, K+ 18)+
+ .125 *V(7,33, K+ 18)+.125 *V(8,33, K+ 18)+
+ .125*V(7,34, K+ 18) +.125 *V(8,34, K+ 18) +
+ .125*V(7,35, K+ 18) +.125*V(8,35,K+ 18)
N E( 1,3, K)=.125 *V(9,32,K+ 18)+ .125*V( 10,32, K+ 18) +
+ .125 *V(9,33, K+ 18)+.125N( 10,33, K+ 18) +
+ .125 *V(9,34, K+ 18) +.125*V(10,34,K+ 18)+
+ .125*V(9,35, K+ 18)+.125*V( 10,35, K+ 18)
SW(1,3, K)=.25*V(7,36,K+ 18) +.25*V(8,36, K+ 18)+
+ .25 *V(7,37 , K+ 18) + .25*V(8,37, K+ 18)
S E(1,3,K) =.25*V(9,36, K+ 18)+.25*V( 10,36, K+ 18) +
+ .25 *V(9,37, K+ 18)+.25*V( 10,37, K+ 18)
NW(2,3,K)=.0625*V(11,32,K+18)+.0625*V(12,32,K+18)+
+ .0625*V(13,32, K+ 18) +.0625*V(14,32, K+ 18) + ,
+ 0625*V(1 1,33, K+ 18)+.0625*V(12,33, K+ 18)+
+ .0625*V(13,33, K+ 18) +.0625*V(14,33, K+ 18) +
+ .0625*V(11,34,K+18)+.0625*V(12,34,K+18)+ ,
+ .0625*V(13,34, K+ 18) +.0625 *V(14,34, K+ 18) +
+ .0625*V( 1 1,35,K+ 18)+.0625 *V(12,35, K+ 18) +
+ .0625*V(13,35, K+ 18) +.0625 *V(14,35,K+ 18)
N E(2,3, K)=.125*V( 15,32, K+ 18)+.125 *V(16,32, K+ 18) +
+ .125*V( 15,33, K+ 18) + .125 *V(16,33, K+ 18)+
+ .125*V(15,34,K+ 18)+.125*V(16,34,K+ 18)+
+ .125*V(15,35,K+ 18)+.125*V(16.35,K+ 18)
SW(2,3, K) =.125*V( 1 1,36, K+ 18) + .125*V( 12,36, K+ 18)+
+ .125*V( 13,36, K+ 18)+.125*V(14,36, K+ 18) +
+ .125 *V( 1 1,37, K+ 18)+ .125*V(12,37, K+ 18) + l
+ .125*V(13,37,K+ 18)+.125*V(14,37,K+ 18)
S E(2,3, K)=.25*V(15,36, K+ 18)+.25*V(16,36,K+ 18)+
+ .25*V(15,37, K+ 18)+.25*V(16,37,K+ 18)
NW(4,3, K)=.0625*V(25,32,K+ 18)+.0625*V(24,32, K+ 18) +
+ .0625*V(25,32,K+18)+.0625*V(26,32,K+18)+
+ .0625*V(23,33, K+ 18)+.0625 *V(24,33, K+ 18)+ l
+ .0625*V(25,33,K+18)+.0625*V(26,33,K+18)+
+ .0625*V(23,34, K+ 18) +.0625*V(24,34, K+ 18) +
+ .0625*V(25,34,K+18)+,0625N(26,34,K+18)+ l
+ .0625*V(23,35, K+ 18)+.0625N(24,35, K+ 18)+
+ .0625*V(25,35,K+ 18)+.0625N(26,35, K+ 18)
N E(4,3,K)=.125*V(27,32,K+ 18)+.125*V(28,32,K+ 18)+
+ .125*V(27,33, K+ 18)+.125N(28,33, K+ 18)+
i
+ .125*V(27,34, K+ 18)+.125*V(28,34, K+ 18)+
l + .125*V(27,35,K+18)+.125*V(28,35,K+18) l SW(4,3,K)=.125N(23,36,K+18)+.125*V(24,36,K+18)+
+ .125*V(25,36,K+ 18)+.125 *V(26,36, K+ 18) +
+ .125*V(23,37, K+ 18)+.125*V(24,37, K+ 18) +
+ .125*V(25,37,K+18)+.125*V(26,37,K+18)
S E(4,3, K)=.25*V(27,36,K+ 18) + .25 *V(28,36, K+ 18)+
+ .25*V(27,37,K+ 18)+.25*V(28,37, K+ 18)
NW(5,3,K)=.0625*V(29,32,K+18)+.0625*V(30,32,K+18)+ ;
+ .0625*V(31,32,K+ 18) +.0625N(32,32, K+ 18)+ l
+ .0625*V(29,33, K+ 18)+.0625*V(30,33, K+ 18) + j j + .0625*V(31,33,K+ 18)+.0625*V(32,33, K+ 18) + !
+ .0625*V(29,34,K+ 18)+.0625*V(30,34,K+ 18) +
+ .0625*V(31,34,K+ 18)+.0625*V(32,34, K+ 18) +
+ .0625N(29,35,K+ 18)+.0625*V(30,35,K+ 18)+
l l
I 4 I 119 l
+ .0625 *V( 31,35, K+ 18)+ .0625 *V(32,35, K+ 18)
N E(5,3,K)=.125*V(33,32,K+ 18)+,125*V(34,32,K+ 18) +
+ .125*V(33,33,K+ 18)+.125*V(34,33,K+ 18) + + .125*V(33,34, K + 18) +.125*V(34,34, K+ 18) + + .125 *V(33,35, K+ 18) +.125*V(34,35, K+ 18)
SW(5,3,K)=.125*V(29,36,K+ 18)+.125*V(30,36,K+ 18)+
+ .125*V(31,36,K+ 18)+.125*V(32,36,K+ 18)+ * + .125*V(29,37,K+18)+.125*V(30,37,K+18)+ + .125*V(31,37, K+ 18)+.125*V(32,37, K+ 18)
S E(5,3, K)= .25*V(33,36, K+ 18) +.25*V(34,36, K+ 18)+
+ .25 *V(33,37, K+ 18) +.25*V(34,37, K+ 10)
NW( 5,2, K) = .125*V(29,38, K+ 18)+.125*V(30,38, K+ 18) +
+ .126*V(31,38, K+ 18) +.125*V(32,38, K+ 18)+ + .125 *V(29,39, K+ 18)+ .125*V(30,39, K+ 18) + +~ .125*V(31,39,K+ 18) +.125*V(32,39, K+ 18)
N E(5,2, K)=.25*V( 33,38, K+ 18) +.25*V(34,38,K+ 18) +
+ .25*V(33,39, K+ 18) +.25*V(34,39, K+ 18)
SW(5,2, K) =.125*V(29,40, K+ 18) +.125*V(30,40, K+ 18)+
+ .125*V(31,40, K+ 18)+.125*V(32,40, K+ 18) + + .125 *V(29,41,K+ 18)+.125*V(30,41, K+ 18)+ + .125 *V(31,41, K+ 18) +.125*V(32,41, K+ 18) ,
S E(5,2,K)=.25*V(33,40,K+ 18)+.25*V(34,40,K+ 18)+ ,
+ .25 *V(33,41, K+ 18) +.25*V(34,41, K+ 18)
ENDDO DO J=7,2,-1 DO i=1,5 APNW(1,J)=0.0 APNE(1,J)=0.0 APSW(1,J)=0.0 APSE (1,J)=0.0 ENDDO ENDDO l DO J=7,2,-1 DO l=1,5 DO K=1,15 APNW(1,J)=APNW(1,J)+NW(1,J,K)/15. AP N E(1,J)=APN E(1,J)+ N E(1,J,K)/15. l APSW(1,J)=APSW(1,J)+SW(1,J,K)/15. APSE (1,J)= APSE (1,J)+SE(1,J,K)/15. l ENDDO ENDDO I ENDDO l DO J=7,2,-1 DO i=1,5 IF(J.EQ.2.AND.I.NE.5)GO TO 10 ) DO K=1,15 J IF(APNW(1,J).NE.O.0)THEN - RPNW(1,J,K)=NW(1,J,K)/APNW(1,J) ELSE RPNW(1,J,K)=0.0 ENDIF ,. IF(APNE(1,J).NE.0.0)THEN RPNE(1,J,K)=NE(1,J,K)/APNE(1,J) l ELSE RPNE(I J,K)=0.0 ENDIF IF(APSW(1,J).NE.0.0)THEN l
i l 120 l i I RPSW(1,J,K)=SW(1,J,K)/APSW(1,J) ) ELSE RPSW(1,J,K)=0.0 ENDIF IF(APSE (1,J).NE.O.0)THEN , RPSE(1,J,K)=SE(1,J,K)/ APSE (1,J) i ELSE l RPSE(1,J,K)=0.0 ENDIF , ENDDO 10 ENDDO I ENDDO PD=0.0 ; I NCOUNT=0 l DO J=7,2,1 DO l=1,5 ) IF(J.EO.2.AND.l.NE.5)GO TO 30 IF(1.EO.3.AND.(J.EQ.7.OR.J.EQ.3))GO TO 30 PD=PD+APNE(1,J)/86. PD=PD+ APSE (1,J)/86. IF(l.EO.1.AND.(J.EO.6.OR.J.EQ.3))GO TO 30 IF(1.EO.5.AND.J.EO.3)GO TO 20 PD=PD+APNW(1,J)/86. 20 IF((l.EO.2.OR.I.EQ.4).AND.(J.EO.6.OR.J.EQ.4))GO TO 30 i IF(1.EQ.3.AND.J.EO.5)GO TO 30 : PD=PD+APSW(1,J)/86. 30 ENDDO ENDDO XCOOR(1)='B' ' XCOOR(2)='C' XCOOR(3)='D' l l XCOOR(4)='E' l XCOOR(5)='F' rrmaxrad=0.0 DO J=7,2,-1 DO l=1,5 IF(1.EO.3.AND.(J.EQ.7.OR.J.EQ.3))GO TO 40 IF(J.EQ.2.AND.I.NE.5)GO TO 40 PDNW=APNW(1,J)/PD if(PDNW.gt.rmaxrad)then rmaxrad=PDNW - maxdir='NW maxcord=xcoor(i) maxpos=j endif PDNE=APNE(1,J)/PD , if(PDNE.gt.rmaxrad)then ' rmaxrad=PDNE l maxdir='NE' l maxcord=xcoor(i) maxpos=j l endif
- PDSW=APSW(1,J)/PD l if(PDSW.gt.rmaxrad)then ,
rmaxrad=PDSW j maxdir='SW
- maxcord=xcoor(i) 1
121 maxpos=j endif PDSE= APSE (1,J)/PD if(PDSE.gt.rmaxrad)then rmaxrad=PDSE maxdir='SE' maxcord=xcoor(i) maxpos=j endif WRITE (20,250)XCOOR(I),J,PDNW,PDNE,PDSW,PDSE ' 40 ENDDO ENDDO write (20,210)maxcord,maxpos,maxdir,rmaxrad DO J=7,2,-1 DO l=1,5 IF(J.EQ.2.AND.l.NE.5)GO TO 50 WRITE (20,300)XCOOR(1),J , DO K=15,1,-1 WRITE (20,400) K, R P NW(1,J, K), R P N E(1,J, K),
+ RPSW(1,J,K),RPSE(1,J,K)
ENDDO 50 ENDDO ENDDO 100 FORMAT (////) 200 FORMAT (10X,9(2X,E11.5)) , 210 format (/,'THE MAXIMUM PEAKING FACTOR OCCURED IN POSITION:',A1,
+11,/,'IN THE ',A2,' DIRECTION WITH A VALUE OF ',F11.8) 250 FORMAT (/,' FUEL PIN POSITION : ',A1,11,/,
+ 'NW : ',F11.8,' NE : ',F11.8,/,
+ 'SW : ',F11.8,' SE : ',F11.8) ,
300 FORMAT (/,' FUEL PIN POSITION : ',A1,li,/,
+ ' Z',7X,'NW,11 X,'N E',11 X,'SW,1 1 X,'S E',/,
+'- ') .
400 FORMAT (12,4(2X,F11.8)) ! C STOP END AXPOW ) l l The second processing code is AXPOW. It reads the output file from Power 2 and extrapolates and interpolates the axial distributions from the nodes generated in DIF3D to the node points necessary for PARET and NCTRIGA. It then prints the output in the form ofinput cards for both PARET and NCTRIGA. The source listing is shown below.
, program axpow J C****************************************************************
. c this program will read in a 15 node axial power distribution c and generate the proper axial distribution data for NCTRIGA and l 4
l 1 1 l 122 l l c PARET. c l l c This was written by Brad Rearden l c for the Nuclear Science Center > c Texas A&M University c September 28,1994 c ; c Modified on October 19,1994
, c This version allows for the input of hot channel l c data and a radial peaking factor to ease the input l c for the new pulsing modelin PARET c
l c Modified on April 25,1995 l c This version will read all of the channels in the core l c and print PARET and NCTRIGA cards for each. This is c intended to eliminate many intermediate steps that c were very time consuming in testing the temperatures e at various core positions c'*"""**"****""****"*""****"***"*"**"**"*"""**** ' c23456789 123456789 123456789 123456789 123456789 123456789 123456789 12 , real f(0:20),xx(0:14),tx(0:20),tf(0:20) real peaks (24,4), shape (0:14,4) character position (24)*2, check *2, direction (4)*2 ; data position / '87', 'C7', 'E7', 'F7', 'B6', 'C6', '
+'D6', 'E6', 'F6', 'B5', 'C5', 'D5', 'E5', 'F5', 'B4',
+'C4', 'D4', 'E4', 'F4', '83', 'C3', 'E3', 'F3', 'F2' / !
data xx / 0.0127,0.0381,0.0635,0.0889,0.1143,
+0.1397, 0.1651, 0.1905, 0.2159, 0.2413, 0.2667,
+0.2921, 0.3175, 0.3429, 0.3683 /
l data direction / 'NW, 'NE', 'SW, 'SE'l l o pe n( u nit = 10, file =' power.dat',statu s='old',fo rm='fo rm atte d') o pe n(u nit = 11, file ='axpow.ou t',statu s='n ew', form ='fo rmatted') c"""""""""""*"""""""*"""*"""""""""" c c This section will read in the radial peaking factors
- c for each fuel rod
! c ! c""*"""*"""""*""*""*"""""""""""""""" lave =1 jave=1 ihot=1 Jhot=1 do i=1,24
)
read (10,8020)(peaks (1,j),j=1,4) do j=1,4 if(abs (peaks (1,j)-1.0).It. abs (peaks (lave.jave)-1.0))then iave=i j jave=j j endif I if(peaks (i,j).gt. peaks (ihot,jhot)) then ihot=i jhot=j endif enddo )
- 8020 format (/,/,6x,f10.8,9x,f10.8,1,6x,f10.8,9x,f10.8) {
enddo 1 I
. write (11,9111) position (lave), direction (jave), peaks (lave,Jave),
l l I 1
123
+ position (ihot), direction (jhot), peaks (ihot.jhot) 9111 format (/,/,'The average rod is ',A,A,' with a value of ',
+ f6.4,/,/,
+'The hot rod is ',A,A,' with a value of ',f6.4) read (10,8021) l 8021 format (2(/)) '
c**"""***""""""""******"********"*"*"******"""" C c Loop over the assemblies and generate the shapes for c each rod c c*"""""""""""""""""""""""""""""""" , do k=1,24 read (10,8030) check 8030 format (/,20x,A) if(check.eq. position (k))then read (10,8040) 8040 format (1(/)) else read (10,8050) 8050 format (20(/)) endif i do !=0,14 l read (10,*)idum,(shape (1,j),j=1,4) enddo ' c""*""""""""""""*****""**"""""""""""*" c l c Loop over the rods in the current assembly to generate j c the power distributions c c*""""*""""***""""*""*"""""""*"*""*"""* } do j=1,4 i do l=0,14 l f(l)= shape (1,j) enddo c*""""*"""""*"""""""""*"""""""""""*" C c Generate the data for use in NCTRIGA and PARET for the ! c average channel c l c"""*"""*"""""""*"""""""""""""""*""* tx(0)=0.0 l tf(0)=f(0)+(f(1)-f(0))/(xx(1)-xx(0))*( 0127) l tx(15)=1.0 tf(15)=f(13)+(f(14)-f(13))/(xx(14)-xx(13))*(.0381) l do 20 i=1,14 i tx(i)=real(i)/15.0 ! call neville(14,tf(i),tx(i)*.381,xx,f) 20 continue c'**"""""""""""**"****""*""""""""*""""" C c Find the average value of the power distibution , c
! c""*"""""""""""""*"""""*"""""""""""
ave =0.0 da 70 i=0,14 l s xl=tx(i)*.381 ;
)
I i l i l l
i ! 124 , xip1=tx(i+1)*0.381 l pown ew= tf(i) + (tf(i+ 1 )-tf(i))/(xip 1 -xi)*((xip 1 -xi)/2.0) i ave = ave +pownew 70 continue I ave = ave /15.0 c print *,' the average is', ave i c"""*"""**""""*"*""""*"*"****"**"***"*""*"* ! C l c Normalize this average to 1.00 ) l c , ! c*""*"""""""*""""""""*"""*"*"""""""*" + l c print *,' normalizing' , l do 100 i=0,15 l tf(i)=tf(i)/ ave ' 100 continue ave =0.0 , do 110 i=0,14
)
xi=tx(i)*.381 xip1=tx(i+1)*0.381 ! pownew=tf(i)+(tf(i+1)-tf(i))/(xip1 xi)*((xip1-xi)/2.0) ave = ave +pownew 110 continue ave = ave /15.0 C***************************************************************** , C c Print the output c , C***************************************************************** write (11,9060) position (k), direction (j) , 9060 format (/,/,/,The following data is formatted for use in',
+ ' NCTRIGA.',/,' This is for ',A,A,/,/) !
write (11,9030)(tx(i),i=0,15) ; 9030 format (2x,'ZR=',20(f6.4,',')) ! write (11,9040)(tf(i),1=0,15) 9040 format (2x,'CVZ=',20(f6.4,','),/,/,/) write (11,9070) position (k), direction (J) 9070 format (/,/,/,
- The following data is formatted for use in PARET,/
+ ' This is for ',A,A,/)
write (11,9050)2,tf(0)* peaks (k,j),1.0.1.0,1.0 do 120 i=1,7 write (11,9050)i+2,f(i)* peaks (k,j),1.0,1.0,1.0 120 continue - do 130 i=8,13 ! write (11,9080)i+2,f(i)* peaks (k,j),1.0,1.0,1.0 ' 130 continue write (11,9080)16,tf(15)* peaks (k,j),1.0.1.0,1.0 9050 format (*510',i1,',',4x,f6.4,6x,f4.2,2(8x,f4.2)) l 9080 format ('51',i2,',',4x,f6.4,6x,f4.2,2(8x,f4.2)) l if((k.eq.iave).and.(J.eq.jave))then write (11,9051)2,tf(0)* peaks (k,j),1.0,1.0,1.0 { do i=1,7 ' write (11,9051)i+2,f(i)* peaks (k,j),1.0,1.0,1.0 enddo do i=8,13
- write (11,9081)i+2,f(i)* peaks (k,j),1.0,1.0,1,0 enddo l
l l I , 125 I l write (11,9081)16,tf(15)* peaks (k.j),1.0,1.0,1.0 9051 format ('520',li ,',',4 x,fG.4,6x,f4. 2,2(8x,f4.2)) 9081 fo rmat('52',12,',',4 x,f6.4,6 x,f4.2,2 (8 x,f4.2)) ; endif enddo I enddo , stop end , l c"""""""""*"**"*"*"*"***"""*"""*""*"*"" - c*""""*""""""""***""*"*""""""""""""" subroutine neville(n,px,x,xx,f) c This subroutine will solve for the value of P(x) c using Neville;s method l real f(0:20),xx(0:20),q(0:20,0:20) , do 10 i=0,n : q(1,0)=f(i) l 10 continue i do 30 i=1,n . do 20 j=1,1 l q(1,j)=((x-xx(i-j))*q(i,j-1)-(x-xx(i))*q(1-1,j-1))/(xx(i) !
-xx(i-j))
20 continue 30 continue px=q(n,n) retum end To verify the accuracy of the output from AXPOW, a plot of the data printed for a fuel channel is shown in Figure 55. This shows the total peaking factors generated for i PARET and NCTRIGA, along with the data from DIF3D. 2.5
+NCTRIGA 2 g 0% PARET I
-dr--DIF3D
[ 1.5 - E i u h 0.5 - 0 0 0.05 0.1 0.15 02 025 0.3 0.35 04 Atlal Position (m) l f l Fig. 55. Comparison of axial distributions generated with AXPOW. i ! j
l l l 126 1 It can be seen from Figure 55 that AXPOW calculates accurately interpolations and extrapolations of the data. I 1 GETFLUX The code GETFLUX reads output files from WIMSD4/m and prints input cards for ! PARET with the within rod radial power shape. The source listing is shown below. program getflux c l c This program will read the output from WIMSD4/m and c calculate the thermal flux for each region throughout i c the cell ! c c This was written in FORTRAN 77 c by Brad Rearden c Nuclear Science Center c Texas A&M University , c June 22.1995 c23456789 123456789 123456789 123456789 123456789 123456789 123456789 12 real fiux(15,23), thflux(15),xnode(12),xcent(11),pflux(8) ! real xparet(8),voicell(8) data xnode /0.0, 0.16, 0.23, 0.62, 0.85,1.02,1.17,
+ 1.31,1.43,1.54,1.64,1.74 /
data voicell / 0.149, .448, .743,1.05,1.34,1.64,1.95,2.24 / character test *39 open(unit =10,fi e='wims.out', status ='old', form =' formatted') o pe n(unit = 11, file =' flux.out',statu s='new', form =' formatted') c initialize the thermal flux to zero do i=1,8 i voicell(i)=voicell(i)/10000.0*0.0254 enddo do k=1,35 do j=1,15 thflux(j)=0.0 enddo icheck=0 i c read in the fluxes dowhile(icheck.eq.0) read (10,9000) test : 9000 format (a) if(test.eq.'1 FLUX DENSITY WITH COLUMN TOTALS ONLY') then icheck=1 , endif enddo read (10,9010) i
127 9010 format (1x,/,1x) do i=1,23 read (10,*) idumb,(flux (j,i),j= 1,8) enddo read (10,9020) 9020 format (1x,/,1x,/,1x,/) do i=1,23 read (10,*) idumb,(flux 0,i),j=9,15) enddo e Find the thermal flux by summing groups 15-23 do j=1,15 do i=15,23 thflux(j)=thflux0)+ flux (j,1) enddo enddo e calculate the location of the center of tha node do i=1,11 xcent(i)=(xnode(i)+xnode(i+ 1))/2.0 enddo c use Neville's algorithm to interpolate the flux to the c node points used in PARET x=0.218/2.0 xparet(1)=x call neville(10,pflux_ node,x,xcent,thflux) pflux(1)=pflux_ node do i=2,8 x=x+0.218 xparet(i)=x call neville(10,pflux_ node,x,xcent,thflux) pflux(i)=pflux_ node enddo c multiply the fluxes by the cell volumes psum=0.0 volsum=0.0 do i=1,8 psum=psum+pflux(i)*voicell(i) enddo pnorm=0.000725926/(psum)*0.032656 volnorm=1.0/(volsum/8) do i=1,8 pflux(i)=pflux(i)*pnorm enddo c Print the thermal flux by region j write (11,9060)k*100 1 9060 format (1x,/,/,/,i6,' MWD bumup') do j=1,8 write (11,9050)],j+1,1,pflux(j) enddo read (10,9000) test enddo i 9050 format ('300',11,', 2.17638-3 ',ii,5x,li,5x,f6.4) stop end C'*******""*********"***"******************"**"******"*** c"""""*"*"""""""**""""""**"""*"""""" subroutine neville(n,px,x,xx,f) c This subroutine will solve for the value of P(x)
128 l c using Neville;s method real f(20),xx(20),q(20,20) do 10 i=1 n+1 q(i,1)=f(i) 10 continue do 30 i=2,n+ 1 do 20 j=2,1 q(1,j)=((x-xx(i-j))*q(i,j-1)-(x-xx(i))*q(1-1,j-1))/(xx(i)
-xx(i-j))
20 continue 30 continue px=q(n+1,n+1) retum end l 1 l
_ ~ . _ _ _ . _ _ _ . _ . . _ _ _ _ . . _ _ . _ _ _ _ . _ _ ._ _.. n a I 129. VITA ' Bradley Thomas Rearden was born in Ft. Worth, Texas on January 23,1971. He is the son of Tom and Sue'Rearden who' currently reside in Granbury, Texas. He attended Arlington High School where he graduated with honors in 1989. He achieved a B.S. in Nuclear Engineering from Texas A&M University in 1993, and 'is currently pursuing a 5 Ph.D in Nuclear Engineering at the same institution. ! During the summers of 1991 and 1992, Brad worked for Texas Utilities in Dallas,
. Texas in the Reactor Engineering group. There, he worked on computer modeling and core following of the Comanche Peak nuclear facility. In the summer of 1993, Brad interned with the Savannah River Site in Aiken, South Carolina where he worked on an input processor for a neutron transport code for the Radiation Shielding group. He began graduate school in August 1993 on a fellowship from the Institute for Nuclear Power Operations. He began working at the Nuclear Science Center in April 1994, and is ,
currently employed there, , His permanent address is the following: , 1006 Indian Creek Ct. , Granbury, TX 76049
)
. _ - - _}}