# Benchmark Analysis on Loss-of-Flow-without-Scram Test of FFTF Using Refined SAC-3D Models

Authors: Siyu Lyu; Daogang Lu; Danting Sui
Journal: Science and Technology of Nuclear Installations
 (2021)

---

## Abstract

The Fast Flux Test Facility (FFTF) is a liquid sodium-cooled nuclear reactor designed by the Westinghouse Electric Corporation for the U.S. Department of Energy. In July 1986, a series of unprotected transients were performed to demonstrate the passive safety of FFTF. Among these, a total of 13 loss-of-flow-without scram (LOFWOS) tests were conducted to confirm the liquid metal reactor safety margins, provide data for computer code validation, and demonstrate the inherent and passive safety benefits of specific design features. In our preliminary work, we have performed relatively coarse modeling of the FFTF. To better predict the transient behavior of FFTF LOFWOS test #13, we modeled it using a more refined thermal-hydraulics model. In this paper, we simulate FFTF LOFWOS test #13 with the system safety analysis code SAC-3D according to the benchmark specifications provided by Argonne National Laboratory (ANL). The simulation range includes the primary and secondary circuits. The reactor core was modeled by the built-in 3D neutronics calculation module and the parallel-channel thermal-hydraulics calculation module. To better predict the reactivity feedback introduced by coolant level variations within the GEMs, a real-time macro cross-section homogenization processing module was developed. The steady-state power distribution was calculated as the transient simulation initial boundary conditions. In general, both the steady-state calculation results and the whole-plant transient behavior predictions are in good agreement with the measured data. The relatively large deviations in transient simulation occur in the outlet temperature predictions of the PIOTA in row 6. It can be preliminarily explained by the reason for neglecting the heat transfer between channels in this model.


---

## Body

## 1. Introduction

In 2017, the International Atomic Energy Agency (IAEA) initialed a Coordinated Research Project (CRP) on benchmark analysis of FFTF loss-of-flow-without scram (LOFWOS) test. A total of 25 organizations from 13 countries participated in this CRP. The overall objective of the CRP is to improve the Member States' analytical capabilities in the field of fast reactor simulation and design. The test defined in this benchmark is LOFWOS test #13, which was initiated at 50% power and 100% flow with the pump pony motors turned off. The conditions of the LOFWOS test, along with the feedback from FFTF's limited free bow core restraint system and the novel passive safety gas expansion modules (GEMs), pose a very challenging and uniquely valuable benchmark exercise.North China Electric Power University (NCEPU) participated in both neutronics and thermal-hydraulics benchmark analysis tasks of the CRP. The fast reactor system analysis code SAC-3D was selected to simulate the steady-state and transient-state for this benchmark. SAC-3D was developed by NCEPU; it was used in fast breeder reactor plant state design and safety analysis [1]. The code was verified using various data experimental measured at mockups with water or sodium coolant and validated using measured data of the CEFR and EBR-II reactors [2]. The code can calculate the process parameters (power, temperature, and flow) of plant transients including normal operating, anticipated operational occurrence (AOO), and design basis accident (DBA). In our preliminary work, we have performed relatively coarse modeling of the FFTF and calculated neutronics and thermal-hydraulics transients [3]. The neutron calculation results are in good agreement, but the thermal-hydraulic transient results deviate significantly. To better predict the transient behavior of FFTF LOFWOS test #13, we modeled it using a more refined thermal-hydraulics model. In this work, we used ERANOS to generate the homogenized macro cross-section library first. Then, the primary and secondary loops were modeled by using SAC-3D. The primary pumps, IHX, and pipes are simulated referring to the technical specifications provided by IAEA. Since the specifications provide the dump heat exchangers (DHXs) sodium outlet temperature and secondary loop mass flow rate as boundary conditions, the pumps and DHXs in the secondary loops were not modeled. Furthermore, to simulate the GEMs properly, we developed a new function module to conduct the macro cross-section homogenization during the transient simulation. The blind calculation phase is now over; most of the measured data has been made available to the benchmark participants. Some data, such as power distribution in a steady state, are still “blind”. For the actual measurement data, we compare the simulation results with the experimental data. For the “blind” part, we compared the simulation results with the results of other participants who applied similar system analysis software for modeling calculations. The neutronics calculation results in both steady-state and transient-state obtained good agreement in general. The trends of the thermal-hydraulics transient calculations are in good agreement, but the values somewhat deviate.


## 2. FFTF LOFWOS Test #13

The Fast Flux Test Facility was a 400 MW-thermal loop type SFR prototype with mixed oxide fuel. The heat was removed from the reactor core by liquid sodium circulating under low pressure. Sodium exited the reactor vessel into one of three primary sodium loops. Intermediate heat exchangers (IHX) separated activated sodium coolant in the primary loops from nonradioactive sodium in the secondary loops. FFTF did not generate electricity, instead of rejecting all heat to the environment via twelve air dump heat exchangers (DHX). Figure1 illustrates the major components in each of the three coolant systems.Figure 1
Overview of FFTF coolant system.LOFWOS test #13 was performed on July 18, 1986. Table1 illustrates the initial conditions of FFTF LOFWOS test #13. It started from 50% power and 100% flow when the three primary pumps were simultaneously tripped. The secondary loop sodium pumps remained operational throughout the test. The plant protection system (PPS) was modified to allow the test to run without the control rods being inserted prematurely.Table 1
Initial conditions of FFTF LOFWOS test #13 [4].
ParameterUnitsValuePowerMWt199.2Core inlet temperatureK590.4Flow through all assemblieskg/s1988.42Shield flow ratekg/s56.32Leakage and bypass flow ratekg/s157.49GEM sodium level1cm221.6Loop 1Loop 2Loop 3Primary pump speedrpm953.1951.5944.4Primary loop flow ratekg/s736.91735.68729.65Secondary pump speedrpm858.2861851.5Secondary loop flow ratekg/s735.18737.59729.47Average DHX Na outlet temperatureK575.8572.6575.4Secondary cold leg temperature2K570.7569.8570.91Elevation relative to the top of inlet nozzle holes.2Initial DHX outlet temperatures were higher than secondary cold leg temperatures, likely due to uncertainties and the location of the thermocouples within the DHX.The FFTF core has 199 hexagonal assemblies. Figure2 shows the assemblies’ arrangement at steady-state prior the test #13, which included eight different types of assemblies: driver fuel assemblies (DFA), in core shim assemblies (ICSAs), reflector assemblies (REFL), control rods (CRs), safety rods (SRs), materials open test assembly (MOTA), fracture mechanics assembly (FMA), and gas expansion modules (GEMs). Among these assemblies, ICSA, MOTA, and FMA are additional types of experimental assemblies, which serve different irradiation and testing purposes. GEMs are nine specially designed assemblies for FFTF. They were introduced to the LOFWOS to create negative reactivity feedback. Each GEM was inserted into an assembly position within the inner row of the radial reflector region. Figure 3 depicts the configuration of GEM. It was sealed at the top and open at the bottom. The GEM utilized its internal 0.0283 m3 volume to trap an argon cover gas bubble when it was first inserted into the sodium [4]. This gas bubble was trapped as long as the assembly remained immersed in sodium. The gas and sodium interface will keep a stable level above the active region's top according to the static sodium head and the pumps' pressure head in normal operating conditions. When the primary pumps are tripped, the pressure head will reduce. Then, the gas bubble in the GEMs expands. As a result, the sodium level drops back down below the active core level, introducing automatic negative reactivity feedback caused by the increasing radial neutron leakage from the core.Figure 2
Core loading for LOFWOS test #13.Figure 3
Gas expansion module assembly [4].To monitor the outlet temperatures of core assemblies during transient conditions, special fast-response thermocouple packages were utilized for a selection of assemblies. These assemblies were referred to as Proximity Instrumented Open Test Assemblies (PIOTAs). Above each PIOTA, a separate duct assembly with multiple thermocouples was installed. The instrument duct was in direct contact with the assembly duct such that sodium discharged from PIOTA assemblies flowed along with thermocouples just above the assembly outlets. Two fast-response PIOTAs were installed for Test #13, and their position is provided in Figure2.


## 3. Calculation Methodology

### 3.1. Neutronics Calculation Model

For neutronics calculation modeling, the core is simulated following the assembly-by-assembly approach. We modeled the FFTF reactor core with the built-in neutronics calculation module of SAC-3D [1, 5]. This module uses the high-order hexagonal nodal expansion method to solve the diffusion equation's formulation in an unlimited number of energy groups [6]. Thus, the self-shielded macroscopic cross-sections are needed as input terms to initiate the calculation. In this benchmark work, the cell/lattice calculation module (named ECCO) integrated with ERANOS [7] was used to generate the macroscopic cross-section library of all materials. As we mentioned before, there are a total of eight types of assemblies within the core. The ICSA, FMA, and MOTA have all different geometries and material specifications, as the ICSA and MOTA are very similar to DFA while, for instance, the FMA is mostly comparable to REFL. Since the output required by the current benchmark analysis does not include the parameters related to the aforementioned assemblies, it was decided to group these three assemblies as DF and REFL.The benchmark specifications provided detailed material compositions of all assemblies; i.e., each assembly had its material composition with 13 axial variations. There are a total of 199 assemblies, which results in a large amount of material requiring cell calculation. To reduce the calculation burden, we grouped all materials into 13 representative mixtures for cell calculations based on the atomic density of the material's significant elements and the assembly type. Table2 shows the cell calculation treatments of all 13 materials. As shown in Table 2, DF3.1 to DF4.2 represent the active zone materials in each of the four types of fuel assemblies. The lower insulator and upper insulator represent the insulator zone materials of driver fuel assemblies separately. The material SR/CR Ab represents the absorber zone of the safety rods and control rods. Sodium CH (SR/CR) and Sodium CH (GEM) represent the materials in the areas of the control rods, safety rods, and GEM assemblies that are filled with sodium coolant, respectively. Similarly, the material GAS represents the material in the gas-filled region of the GEM assemblies. Since the geometry and elemental composition of the reflector assemblies are almost constant within the neutronics modeled range, we consider the reflector assembly as one material.Table 2
Cell calculation treatment of FFTF LOFWOS test #13.
Material NO.XS nameAssembly typeGeometry treatmentXS treatment1DF3.1Driver fuelHeterogeneousCritical2DF4.1Driver fuelHeterogeneousCritical3DF3.2Driver fuelHeterogeneousCritical4DF4.2Driver fuelHeterogeneousCritical5Lower insulatorDriver fuelHomogeneousSubcritical6Upper insulatorDriver fuelHomogeneousSubcritical7SR/CR AbSafety/Control rodHeterogeneousSubcritical8Sodium CH (SR/CR)Safety/Control rodHomogeneousSubcritical9Below AbSafety/Control rodHomogeneousSubcritical10Ref 7/7.2ReflectorHomogeneousSubcritical11Ref 9ReflectorHomogeneousSubcritical12Sodium CH (GEM)GEMHomogeneousSubcritical13GASGEMHomogeneousSubcriticalCell calculations have been performed adopting a 33-energy-group structure and the P1 approximation for scattering treatment and evaluating the thermal expansion and Doppler effect based on the ENDF/B-VII.1.1968-group nuclear data library. For the fissile cells, the heterogeneous hexagonal geometry has been adopted in a 6-step calculation procedure. The heterogeneous hexagonal geometry description has also been used for control and safety rods spatial self-shielded cross section (4-step calculation). The rest nonburnable regions (insulator, reflectors, etc.) foresee a homogeneous geometry description of the cell and use the 2-step calculation procedure. Since the subcritical cells do not have the inner neutron source, we calculate a relative region's spectrum for each of them as the external neutron source. The leakage effects of these subcritical cells are also being considered. Users can define the geometrical buckling of the subcritical cells in the input card of ECCO to take the leakage effects into account. The leakage effects are calculated by the following equations according to the coding manual.(1)B2=58πLtyp2,where B is the buckling value in the correction method Σt+DB2, and Ltyp is the typical dimension of the subcritical region.To update the cross section more accurately according to the material temperature in transient calculations, we performed cell calculations at 500 K, 750 K, and 1000 K, respectively, and obtained two sets of homogeneous macro cross-section libraries. The other material cross-section libraries under temperatures between 500 K and 1000 K will be derived by linear interpolation.The neutron calculation model's computational region covers each assembly in the radial direction. The axial computational region is the area from the bottom to the top of the active region. Each assembly was set as one calculation node in the radial direction. Figure4 is the radial node division schematic diagram of the FFTF reactor core. In the axial direction, the neutronics calculation module of SAC-3D allows users to define the calculation nodes according to their actual needs. As Figure 5 shows, we divided each assembly model into 38 calculation blocks in the present work. Each calculation block is 2.5 cm tall, while 38 blocks have a total of 95 cm height, which respects the active region's actual height. The relative axial positions of the control rod materials with respect to the fuel and GEM assemblies are also shown in Figure 5. By dividing the reactor core calculation model into several calculation nodes in the radial and axial directions, we have 7562 calculation nodes in the neutronics calculation model. Zero flux boundary conditions are used to close the neutronic problem at the core edges.Figure 4
Loading pattern of the calculation model of FFTF in the radial direction.Figure 5
Schematic of calculation nodes division of FFTF in the axial direction.As mentioned earlier, the sodium level in GEMs is higher than the top of the core active region. When the test began, the sodium level in GEMs will drop rapidly to the inlet nozzle region, which is lower than the bottom of the active region. To predict the sodium level, the following correlation was suggested by the technical specifications:(2)L=265−5395042440.13+F2,where F is the core flow rate in percent, and L is the level of sodium above the top of the inlet flow holes.In the transient calculation, sodium and gas within GEMs might be present simultaneously in the same calculation blocks where the gas-liquid interface is located. Therefore, the original cross-sections employed in the static calculations of these blocks need to be specially processed. We solve this problem by conducting the cross-section homogenization during the online feedback reactivity calculation. The homogenization method was shown in Figure6. At first, we assumed that each assembly's boundary condition in the radial direction is a total reflection and kept the axial direction boundary conditions the same as in the core neutronics calculation. Then, we get the neutron flux distribution in the axial direction by conducting the one-dimensional SN (Discrete ordinates method) calculation and do the needed cross-section homogenization as follows [8]:(3)∑i=∫0h′∑i1φdh+∫h′h∑i2φdh∫0hφdh,where Σi are the cross sections consisting of removal, fission, scattering, absorption, and total cross sections, Σi1 and Σi2 are the cross sections from the two neighbor nodes, respectively, h is the height of each node, h′ is the height of some point of the calculation block, and φ is the neutron flux.Figure 6
Cross-section homogenization process.


### 3.2. Thermal-Hydraulics Calculation Model

The design of the reactor was a three-loop primary system where, under low pressure, the three primary loops transport the heat generated within the reactor core to the intermediate heat exchangers and back to the reactor vessel. Each of the three primary loops contains the same components. Given the intrinsic symmetry characteristics of the thermal-hydraulic conditions of the primary loops within the framework of the present benchmark activities, it was decided for the sake of simplicity to implement a single equivalent primary and secondary loop model. For the IHX, primary pump, and pipe models in this work, we took into account their heat transfer area, mass flow rate, and other factors to ensure that they could be equivalent to the IHXs, primary pumps, and pipes in the overall three primary loops. Therefore, the SAC-3D model of the IHXs and primary pumps only consists of a single IHX and a single pump, and a single average pipe layout is defined for the primary loops. As for the primary loop, a single equivalent secondary loop is modeled in SAC-3D.Figure7 shows the layout of the FFTF reactor system simulation model. The primary loop simulation model was divided into two parts: the reactor vessel and the primary loop, which include the hot leg piping, the primary pump, the primary side of IHXs, and the cold leg piping. Hot leg piping runs from the outlet plenum of the reactor vessel to the primary pump and into the IHX primary side. Cold leg piping returns sodium from the IHX to the inlet plenum of the reactor vessel, discharging into the inlet plenum of the reactor vessel. With regard to the secondary loops, only the secondary side of the IHX was modeled, imposing proper boundary conditions (coolant flow rate and temperature) at the secondary hot leg and the secondary cold leg.Figure 7
Simplified layout of FFTF Simulation Model.SAC-3D uses the parallel-channel approach to conduct the reactor vessel thermal-hydraulics calculation, which does not consider the coolant crossflow between the internal components. Simultaneously, the model does not consider the fuel component's axial heat conduction [9]. In the preliminary work, we used only 9 channels to simulate the core flow path and did not consider the bypass flow [3]. This approximate treatment brings large deviations in the thermal-hydraulic transient results. To better correct the bias due to the roughness of the thermal-hydraulic model, a more refined modeling approach was used. In this paper, the reactor vessel model consists of 13 parallel flow channels (11 within the core model) connecting the inlet and outlet plenums. The inlet and outlet plenums are modeled with a zero-dimensional volume module with respect to their effective height and volume. Within the outlet plenum model, the sodium level is set at its nominal value and the argon gas pressure boundary condition is set at the top of the volume. As shown in Figure 8, the reactor core model consists of 11 parallel flow channels; all assemblies are grouped by their coolant flow rate, power distribution, and assembly type. Channels 1 to 6 with channels 10 and 11 represent driver fuel assemblies flow channels, where channels 10 and 11 represent the individual PIOTA channels located in rows 2 and 6, respectively. Figure 9 shows the nodalization scheme of the driver fuel channels. Each of these channels connects to the outlet model and the inlet model from their top and bottom. These models are set to take the loss of pressure caused by the orifice region and wire-wrapped fuel pin bundle into account. As shown in Figure 9, each fuel channel is divided into five calculation points in radial directions. In the axial direction, the number of control volumes is set as the same as the neutronics model to implement neutronics thermal-hydraulics coupling during the transient simulation. Table 3 shows the initial condition comparison between the measurements and SAC-3D model results. Channels 7 represent control/safety rod flow channels. Channels 8 and 9 represent two types of reflector assemblies (REF 7/72/8A are located in the core basket and fed by sodium from the core basket. REF 8 B/9 are located in the support skirt and fed by sodium from the annular plenum), respectively. In addition, several other flow paths draw sodium from the inlet plenum and discharge sodium into the outlet plenum. Leakage and bypass flow paths provide sodium flow to the in-vessel storage region, the gap between the reactor vessel and vessel thermal liner. We simplify these flow paths into two flow channels: In-Vessel Storage Channel and Bypass Channel.Figure 8
SAC-3D channel layout for FFTF core.Figure 9
Thermal-hydraulics model of driver fuel channels.Table 3
Initial core flow distribution comparison.
ChannelAssembly typeMeasured mass flow rate (kg/s)FFTF model Mass flow rate (kg/s)Discrepancy (%)1DF 3.2633.05632.82−0.042DF 4.2152.33152.3−0.023DF 3.1499.07498.54−0.114DF 3.1289.73289.56−0.065DF 4.199.8499.960.126DF 4.1152.33153.230.597SR/CR40.9140.76−0.378REF 7/72/8A52.3852.430.109REF 8 B/923.4823.530.2110PIOTA 224.8224.940.4811PIOTA 620.520.37−0.63In this work, the primary pump is modeled with a pump module based on the pump's characteristics and their nominal operating state. Initial pump speed is adjusted to the benchmark data. The pump speed history is given as input data during the transient. The benchmark specifications provide the pump performance curves for the primary pumps. However, the data for these curves were not available. To predict the pressure head during the transient, the following set of equations of the modified homologous pump theory of Wylie and Streeter were coded into the pump model. These modified correlations provide a good approximation of the dimensionless pump head (H¯).(4)H¯=N¯2+Q¯2Whx,where N¯ is the dimensionless pump speed. Q¯ is the relative volumetric flow rate. x=π+arctanQ¯/N¯,(5)Whx=∑i=06ah,ixi,whose coefficients are provided in Table 4.Table 4
Coefficients for the pump curve.
iah0431.96691−576.614382301.000293−75.4658648.675505−0.260626−0.01596The IHXs were vertically mounted counterflow shell and tube designs. Primary sodium entered the IHX from the outlet of the primary loop hot leg piping and flowed down through the active heat transfer region along the outside of the IHX tubes. Secondary sodium entered the IHX at the top and flowed down through a central downcomer into the lower hemispherical plenum. From this plenum, secondary sodium entered the IHX tubes at the bottom and flowed upward along the inner side of the IHX tubes out of the IHX. In the IHX model, all heat transfer tubes have been simplified as one representative tube. In an attempt to consistently simplify the model, there are several assumptions: (1) there is an ideal mixing of coolant at the inlet and outlet; (2) fully developed convective heat transfer is assumed; (3) no metallic structure of the IHX other than IHX tubes and shell wall is considered; hence, perfect insulation is postulated between primary and secondary sodium in the central downcomer pipe; (4) no heat losses are accounted for through the external IHX vessel to the surrounding environment. As shown in Figure10, there are four radial calculation nodes: secondary coolant, tube metal, primary coolant, and shell wall. The temperature calculation nodes of the heat tube metal and shell wall are defined in the center of the control volume. The temperature calculation nodes of the coolant are defined at the boundary of the control body [9].Figure 10
IHX calculation model.The energy equation of the primary side coolant control volume can be written as(6)ρVPdep,i+1dt=Wpep,i−ep,i+1−HptAptTp,i,i+1−Tt,i−HpshApshTp,i,i+1−Tsh,i.The energy equation of the secondary side coolant control volume can be written as(7)ρVPdes,idt=Wses,i+1−es,i+HstAstTt,i−Ts,i,i+1.The energy equation of the IHX tube metal control volume can be written as(8)Mtct,idTt,idt=HptAptTp,i,i+1−Tt,i−HpstApstTt,i−Ts,i,i+1.The energy equation of the IHX shell wall control volume can be written as(9)Mshcsh,idTsh,idt=HpshApshTp,i,i+1−Tsh,i,where p represents the coolant of the primary side, s represents the coolant of the secondary side, t represents the heat exchanger tube of IHX, sh represents the shell wall of IHX, e is enthalpy, Wis mass flow rate, T is temperature, M is mass, c is specific heat capacity, V is the volume of the control volume, H is the heat transfer coefficient between fluid and structural, A is the heat transfer area between the fluid and the structure, and Ti,i+1 is the average fluid temperature within two neighbor control volumes; it is defined as(10)Ti,i+1=Ti+Ti,i+12.As shown in Figure7, the secondary loop was modeled as two boundary conditions. The DHX outlet coolant temperature values and the secondary loop coolant mass flow rates are imposed as the boundary conditions for the inlet and outlet of the IHX secondary side.


## 3.1. Neutronics Calculation Model

For neutronics calculation modeling, the core is simulated following the assembly-by-assembly approach. We modeled the FFTF reactor core with the built-in neutronics calculation module of SAC-3D [1, 5]. This module uses the high-order hexagonal nodal expansion method to solve the diffusion equation's formulation in an unlimited number of energy groups [6]. Thus, the self-shielded macroscopic cross-sections are needed as input terms to initiate the calculation. In this benchmark work, the cell/lattice calculation module (named ECCO) integrated with ERANOS [7] was used to generate the macroscopic cross-section library of all materials. As we mentioned before, there are a total of eight types of assemblies within the core. The ICSA, FMA, and MOTA have all different geometries and material specifications, as the ICSA and MOTA are very similar to DFA while, for instance, the FMA is mostly comparable to REFL. Since the output required by the current benchmark analysis does not include the parameters related to the aforementioned assemblies, it was decided to group these three assemblies as DF and REFL.The benchmark specifications provided detailed material compositions of all assemblies; i.e., each assembly had its material composition with 13 axial variations. There are a total of 199 assemblies, which results in a large amount of material requiring cell calculation. To reduce the calculation burden, we grouped all materials into 13 representative mixtures for cell calculations based on the atomic density of the material's significant elements and the assembly type. Table2 shows the cell calculation treatments of all 13 materials. As shown in Table 2, DF3.1 to DF4.2 represent the active zone materials in each of the four types of fuel assemblies. The lower insulator and upper insulator represent the insulator zone materials of driver fuel assemblies separately. The material SR/CR Ab represents the absorber zone of the safety rods and control rods. Sodium CH (SR/CR) and Sodium CH (GEM) represent the materials in the areas of the control rods, safety rods, and GEM assemblies that are filled with sodium coolant, respectively. Similarly, the material GAS represents the material in the gas-filled region of the GEM assemblies. Since the geometry and elemental composition of the reflector assemblies are almost constant within the neutronics modeled range, we consider the reflector assembly as one material.Table 2
Cell calculation treatment of FFTF LOFWOS test #13.
Material NO.XS nameAssembly typeGeometry treatmentXS treatment1DF3.1Driver fuelHeterogeneousCritical2DF4.1Driver fuelHeterogeneousCritical3DF3.2Driver fuelHeterogeneousCritical4DF4.2Driver fuelHeterogeneousCritical5Lower insulatorDriver fuelHomogeneousSubcritical6Upper insulatorDriver fuelHomogeneousSubcritical7SR/CR AbSafety/Control rodHeterogeneousSubcritical8Sodium CH (SR/CR)Safety/Control rodHomogeneousSubcritical9Below AbSafety/Control rodHomogeneousSubcritical10Ref 7/7.2ReflectorHomogeneousSubcritical11Ref 9ReflectorHomogeneousSubcritical12Sodium CH (GEM)GEMHomogeneousSubcritical13GASGEMHomogeneousSubcriticalCell calculations have been performed adopting a 33-energy-group structure and the P1 approximation for scattering treatment and evaluating the thermal expansion and Doppler effect based on the ENDF/B-VII.1.1968-group nuclear data library. For the fissile cells, the heterogeneous hexagonal geometry has been adopted in a 6-step calculation procedure. The heterogeneous hexagonal geometry description has also been used for control and safety rods spatial self-shielded cross section (4-step calculation). The rest nonburnable regions (insulator, reflectors, etc.) foresee a homogeneous geometry description of the cell and use the 2-step calculation procedure. Since the subcritical cells do not have the inner neutron source, we calculate a relative region's spectrum for each of them as the external neutron source. The leakage effects of these subcritical cells are also being considered. Users can define the geometrical buckling of the subcritical cells in the input card of ECCO to take the leakage effects into account. The leakage effects are calculated by the following equations according to the coding manual.(1)B2=58πLtyp2,where B is the buckling value in the correction method Σt+DB2, and Ltyp is the typical dimension of the subcritical region.To update the cross section more accurately according to the material temperature in transient calculations, we performed cell calculations at 500 K, 750 K, and 1000 K, respectively, and obtained two sets of homogeneous macro cross-section libraries. The other material cross-section libraries under temperatures between 500 K and 1000 K will be derived by linear interpolation.The neutron calculation model's computational region covers each assembly in the radial direction. The axial computational region is the area from the bottom to the top of the active region. Each assembly was set as one calculation node in the radial direction. Figure4 is the radial node division schematic diagram of the FFTF reactor core. In the axial direction, the neutronics calculation module of SAC-3D allows users to define the calculation nodes according to their actual needs. As Figure 5 shows, we divided each assembly model into 38 calculation blocks in the present work. Each calculation block is 2.5 cm tall, while 38 blocks have a total of 95 cm height, which respects the active region's actual height. The relative axial positions of the control rod materials with respect to the fuel and GEM assemblies are also shown in Figure 5. By dividing the reactor core calculation model into several calculation nodes in the radial and axial directions, we have 7562 calculation nodes in the neutronics calculation model. Zero flux boundary conditions are used to close the neutronic problem at the core edges.Figure 4
Loading pattern of the calculation model of FFTF in the radial direction.Figure 5
Schematic of calculation nodes division of FFTF in the axial direction.As mentioned earlier, the sodium level in GEMs is higher than the top of the core active region. When the test began, the sodium level in GEMs will drop rapidly to the inlet nozzle region, which is lower than the bottom of the active region. To predict the sodium level, the following correlation was suggested by the technical specifications:(2)L=265−5395042440.13+F2,where F is the core flow rate in percent, and L is the level of sodium above the top of the inlet flow holes.In the transient calculation, sodium and gas within GEMs might be present simultaneously in the same calculation blocks where the gas-liquid interface is located. Therefore, the original cross-sections employed in the static calculations of these blocks need to be specially processed. We solve this problem by conducting the cross-section homogenization during the online feedback reactivity calculation. The homogenization method was shown in Figure6. At first, we assumed that each assembly's boundary condition in the radial direction is a total reflection and kept the axial direction boundary conditions the same as in the core neutronics calculation. Then, we get the neutron flux distribution in the axial direction by conducting the one-dimensional SN (Discrete ordinates method) calculation and do the needed cross-section homogenization as follows [8]:(3)∑i=∫0h′∑i1φdh+∫h′h∑i2φdh∫0hφdh,where Σi are the cross sections consisting of removal, fission, scattering, absorption, and total cross sections, Σi1 and Σi2 are the cross sections from the two neighbor nodes, respectively, h is the height of each node, h′ is the height of some point of the calculation block, and φ is the neutron flux.Figure 6
Cross-section homogenization process.


## 3.2. Thermal-Hydraulics Calculation Model

The design of the reactor was a three-loop primary system where, under low pressure, the three primary loops transport the heat generated within the reactor core to the intermediate heat exchangers and back to the reactor vessel. Each of the three primary loops contains the same components. Given the intrinsic symmetry characteristics of the thermal-hydraulic conditions of the primary loops within the framework of the present benchmark activities, it was decided for the sake of simplicity to implement a single equivalent primary and secondary loop model. For the IHX, primary pump, and pipe models in this work, we took into account their heat transfer area, mass flow rate, and other factors to ensure that they could be equivalent to the IHXs, primary pumps, and pipes in the overall three primary loops. Therefore, the SAC-3D model of the IHXs and primary pumps only consists of a single IHX and a single pump, and a single average pipe layout is defined for the primary loops. As for the primary loop, a single equivalent secondary loop is modeled in SAC-3D.Figure7 shows the layout of the FFTF reactor system simulation model. The primary loop simulation model was divided into two parts: the reactor vessel and the primary loop, which include the hot leg piping, the primary pump, the primary side of IHXs, and the cold leg piping. Hot leg piping runs from the outlet plenum of the reactor vessel to the primary pump and into the IHX primary side. Cold leg piping returns sodium from the IHX to the inlet plenum of the reactor vessel, discharging into the inlet plenum of the reactor vessel. With regard to the secondary loops, only the secondary side of the IHX was modeled, imposing proper boundary conditions (coolant flow rate and temperature) at the secondary hot leg and the secondary cold leg.Figure 7
Simplified layout of FFTF Simulation Model.SAC-3D uses the parallel-channel approach to conduct the reactor vessel thermal-hydraulics calculation, which does not consider the coolant crossflow between the internal components. Simultaneously, the model does not consider the fuel component's axial heat conduction [9]. In the preliminary work, we used only 9 channels to simulate the core flow path and did not consider the bypass flow [3]. This approximate treatment brings large deviations in the thermal-hydraulic transient results. To better correct the bias due to the roughness of the thermal-hydraulic model, a more refined modeling approach was used. In this paper, the reactor vessel model consists of 13 parallel flow channels (11 within the core model) connecting the inlet and outlet plenums. The inlet and outlet plenums are modeled with a zero-dimensional volume module with respect to their effective height and volume. Within the outlet plenum model, the sodium level is set at its nominal value and the argon gas pressure boundary condition is set at the top of the volume. As shown in Figure 8, the reactor core model consists of 11 parallel flow channels; all assemblies are grouped by their coolant flow rate, power distribution, and assembly type. Channels 1 to 6 with channels 10 and 11 represent driver fuel assemblies flow channels, where channels 10 and 11 represent the individual PIOTA channels located in rows 2 and 6, respectively. Figure 9 shows the nodalization scheme of the driver fuel channels. Each of these channels connects to the outlet model and the inlet model from their top and bottom. These models are set to take the loss of pressure caused by the orifice region and wire-wrapped fuel pin bundle into account. As shown in Figure 9, each fuel channel is divided into five calculation points in radial directions. In the axial direction, the number of control volumes is set as the same as the neutronics model to implement neutronics thermal-hydraulics coupling during the transient simulation. Table 3 shows the initial condition comparison between the measurements and SAC-3D model results. Channels 7 represent control/safety rod flow channels. Channels 8 and 9 represent two types of reflector assemblies (REF 7/72/8A are located in the core basket and fed by sodium from the core basket. REF 8 B/9 are located in the support skirt and fed by sodium from the annular plenum), respectively. In addition, several other flow paths draw sodium from the inlet plenum and discharge sodium into the outlet plenum. Leakage and bypass flow paths provide sodium flow to the in-vessel storage region, the gap between the reactor vessel and vessel thermal liner. We simplify these flow paths into two flow channels: In-Vessel Storage Channel and Bypass Channel.Figure 8
SAC-3D channel layout for FFTF core.Figure 9
Thermal-hydraulics model of driver fuel channels.Table 3
Initial core flow distribution comparison.
ChannelAssembly typeMeasured mass flow rate (kg/s)FFTF model Mass flow rate (kg/s)Discrepancy (%)1DF 3.2633.05632.82−0.042DF 4.2152.33152.3−0.023DF 3.1499.07498.54−0.114DF 3.1289.73289.56−0.065DF 4.199.8499.960.126DF 4.1152.33153.230.597SR/CR40.9140.76−0.378REF 7/72/8A52.3852.430.109REF 8 B/923.4823.530.2110PIOTA 224.8224.940.4811PIOTA 620.520.37−0.63In this work, the primary pump is modeled with a pump module based on the pump's characteristics and their nominal operating state. Initial pump speed is adjusted to the benchmark data. The pump speed history is given as input data during the transient. The benchmark specifications provide the pump performance curves for the primary pumps. However, the data for these curves were not available. To predict the pressure head during the transient, the following set of equations of the modified homologous pump theory of Wylie and Streeter were coded into the pump model. These modified correlations provide a good approximation of the dimensionless pump head (H¯).(4)H¯=N¯2+Q¯2Whx,where N¯ is the dimensionless pump speed. Q¯ is the relative volumetric flow rate. x=π+arctanQ¯/N¯,(5)Whx=∑i=06ah,ixi,whose coefficients are provided in Table 4.Table 4
Coefficients for the pump curve.
iah0431.96691−576.614382301.000293−75.4658648.675505−0.260626−0.01596The IHXs were vertically mounted counterflow shell and tube designs. Primary sodium entered the IHX from the outlet of the primary loop hot leg piping and flowed down through the active heat transfer region along the outside of the IHX tubes. Secondary sodium entered the IHX at the top and flowed down through a central downcomer into the lower hemispherical plenum. From this plenum, secondary sodium entered the IHX tubes at the bottom and flowed upward along the inner side of the IHX tubes out of the IHX. In the IHX model, all heat transfer tubes have been simplified as one representative tube. In an attempt to consistently simplify the model, there are several assumptions: (1) there is an ideal mixing of coolant at the inlet and outlet; (2) fully developed convective heat transfer is assumed; (3) no metallic structure of the IHX other than IHX tubes and shell wall is considered; hence, perfect insulation is postulated between primary and secondary sodium in the central downcomer pipe; (4) no heat losses are accounted for through the external IHX vessel to the surrounding environment. As shown in Figure10, there are four radial calculation nodes: secondary coolant, tube metal, primary coolant, and shell wall. The temperature calculation nodes of the heat tube metal and shell wall are defined in the center of the control volume. The temperature calculation nodes of the coolant are defined at the boundary of the control body [9].Figure 10
IHX calculation model.The energy equation of the primary side coolant control volume can be written as(6)ρVPdep,i+1dt=Wpep,i−ep,i+1−HptAptTp,i,i+1−Tt,i−HpshApshTp,i,i+1−Tsh,i.The energy equation of the secondary side coolant control volume can be written as(7)ρVPdes,idt=Wses,i+1−es,i+HstAstTt,i−Ts,i,i+1.The energy equation of the IHX tube metal control volume can be written as(8)Mtct,idTt,idt=HptAptTp,i,i+1−Tt,i−HpstApstTt,i−Ts,i,i+1.The energy equation of the IHX shell wall control volume can be written as(9)Mshcsh,idTsh,idt=HpshApshTp,i,i+1−Tsh,i,where p represents the coolant of the primary side, s represents the coolant of the secondary side, t represents the heat exchanger tube of IHX, sh represents the shell wall of IHX, e is enthalpy, Wis mass flow rate, T is temperature, M is mass, c is specific heat capacity, V is the volume of the control volume, H is the heat transfer coefficient between fluid and structural, A is the heat transfer area between the fluid and the structure, and Ti,i+1 is the average fluid temperature within two neighbor control volumes; it is defined as(10)Ti,i+1=Ti+Ti,i+12.As shown in Figure7, the secondary loop was modeled as two boundary conditions. The DHX outlet coolant temperature values and the secondary loop coolant mass flow rates are imposed as the boundary conditions for the inlet and outlet of the IHX secondary side.


## 4. Results and Analysis

Table5 summarizes the comparison of measurements and calculation results of the initial plant conditions. Reactor power was slightly less than 50% and the core inlet temperature was reduced by nearly 80℉from the nominal core inlet temperature. This was done to ensure that the peak fuel temperatures did not rise too high to limit the impact of the test on the durability of the fuel assemblies.Table 5
Comparison of benchmark and SAC-3D steady-state parameters.
ParameterUnitMeasured valueCalculated valueDiscrepancy (%)Core powerMW199.2199.20.00Core inlet temperatureK590.4590.870.08Total core flow ratekg/s2202.242197.98−0.19Mass flow rate through all assemblieskg/s1988.441989.340.05GEM sodium level1cm221.6221.60.00Primary hot leg temperatureK—678.18—Primary cold leg temperatureK—605.89—Secondary hot leg temperatureK—656.56—Secondary cold leg temperatureK583.06583.780.121Elevation relative to the top of inlet nozzle holes.Note that, to date, the IAEA's FFTF CRP project has just completed the blind phase, and ANL has released experimental measurements of the part of the calculated FFTF reactor parameters. Therefore, the following results partly compare the simulated calculations with the real measurements and partly compare with other participants' calculations. Because different participants use different calculation methods and models, we only make data deviation comparisons without deviation analysis for the results without measurement data.Table6 lists the main neutron results of FFTF LOFWOS test #13 of 10 participants. Most of the results, like the neutron multiplication factor, delayed neutron fraction, Doppler coefficient, fuel density coefficient, etc., have a relatively good agreement among participants. After removing results that deviated from the mean by more than two times the standard deviation, the variation coefficient (δ/Avg) was distributed between 1% and 27%. Significant deviations appear in the structure density and sodium density coefficients' results, which have a 542% and 98% variation coefficient, respectively.Table 6
Main results of FFTF LOFWOS test #13 neutronics benchmark analysis
ArgonneIGCARINESTIPPEJAEAKITKTHPSIROMENCEPUNeutron multiplication factor1.0000.9980.9990.9921.0170.9981.0221.0061.0000.998Delayed neutron fraction0.0030.0030.0070.0030.0030.0040.0030.0030.0030.003Axial expansion coefficient (pcm/°C)−0.32−0.23——−0.32−−0.3−0.22−0.48−0.3Radial expansion coefficient (pcm/°C)−1−1.22——−1−−0.93−1.52−5.87−0.94Fuel Doppler constant (pcm)−629−508——−634−509−564−658−688−524Fuel density coefficient (pcm)−1.36−1.45——−0.45—−1.36−1.36−1.4−1.32Structure density coefficient (pcm)−0.120.2——0.03—0.10.04−0.1−0.01Sodium density coefficient (pcm/°C)−0.35−0.91——−0.410.09−0.94−0.27−1.91−0.37Control and safety rods (pcm)−11849——−939610800—−11540−11823−12773−8343Gas expansion modules (pcm)−442−498—−516−489−448420−475−1201−782The power distribution of steady-state estimated by SAC-3D is shown in Figure11. The axially integrated power of each assembly ranges from 1.55 MW to 3.75 MW. Figure 12 illustrates the variation coefficient map obtained by comparing the integrated power simulation results with JAEA and ANL's results. The largest variation coefficient is 3.27%, which occurred in a fuel assembly type 4.2D. Besides, we note that other assemblies with large coefficients of variation are located around this fuel assembly. Since we do not have the power distribution's measured data, we cannot objectively assess whether these biases are significant. By comparing the codes applied by JAEA, ANL, and our side, we believe that these biases are due to, at least part of, the modeling methods differences and different cross-section libraries.Figure 11
Power distribution map of FFTF core.Figure 12
Power distribution calculation variation coefficient map of FFTF core.Figures13–15 depict the primary loop's mass flow rate's transient profile. It can be seen that the SAC-3D very well predicts the pump coast down. Minor differences exist during the natural circulation where the mass flow rates are slightly lower than the experimental data. As shown in Figure 16, this is due to the underestimation of total power.Figure 13
Mass flow rate in loop #1.Figure 14
Mass flow rate in loop #1 (Low range).Figure 15
Mass flow rate in loop #1 (0 s–100 s).Figure 16
Reactor power of FFTF LOFWOS test #13 (low range).Figures16–18 show the transient reactor power profiles. Figures 19 and 20 are the transient reactivity profile. At 0 s, the primary pump is tripped, and the core flow rate decreases rapidly. Due to the rapid decrease of the core cycling capacity, the assemblies' temperature starts to rise and introduce negative net reactivity to bring down the reactor power. Starting from 5 s, the coolant level in the GEMs drops rapidly and the neutron leakage rate increases significantly, introducing a remarkable negative reactivity. The fission power decreases rapidly and the total power drops to 42% of the initial power at 10 s. At 20 s, the introduced negative reactivity increases at a slower rate as the core flow rate and the coolant level within the GEMs begin to decrease at a slower rate. As shown in Figure 14, at approximately 90 s, the natural cycle is established. The coolant level within GEMs begins to stabilize and the negative reactivity it introduces also stabilizes. Figure 16 shows that the fission power is overestimated, due to the underestimation of the negative net reactivity. In the long term, the fission power is close to zero, and the total power is contributed mainly by the decay power. The underestimation of the decay power leads to an underestimation of the total power after 200 seconds. As shown in Figure 19, although the fuel temperature decreases introduced by the positive Doppler's reactivity feedback and positive axial expansion reactivity feedback, the net reactivity still maintains a significant negative value due to the significant negative reactivity introduced by GEMs.Figure 17
Reactor power of FFTF LOFWOS test #13.Figure 18
Reactor power of FFTF LOFWOS test #13 (0 s–40 s).Figure 19
Reactivity transient profile of FFTF LOFWOS test #13.Figure 20
Reactivity Transient Profile of FFTF LOFWOS test #13 (0 s–50 s).Figures21 and 22 display the transient outlet temperature profiles of PIOTAs. The SAC-3D results show a very similar trend compared to the experimental data, with two peaks occurring during the entire duration of the transient. In the beginning, the outlet temperature rises rapidly as the core flow rate drops rapidly as shown in Figure 13. The power decreases slowly because of the negative reactivity feedback introduced by the control rod driveline expansion and radial expansion. The power-to-flow ratio reaches the first peak value as shown in Figure 23. Thus, the PIOTAs' outlet temperatures increase rapidly and reach the first peak value. After that, the rapid drop of the coolant level within the GEMs introduced a large amount of negative reactivity, causing the total power to drop rapidly. Meanwhile, the decreasing rate of mass flow rate starts to slow down. The PIOTAs' outlet temperatures start to decrease because the power-to-flow ratio rapidly decreases. Starting from about 20 s, the decrease in total power starts to slow down due to the rebound of the negative net reactivity as shown in Figure 20. Therefore, the power-to-flow ratio increases which make the PIOTAs’ temperatures rise again. At about 110 s, the PIOTAs’ temperatures reach the second peak value after the power-to-flow ratio reaches the second peak value. Since then, the natural circulation starts to be established and maintained at about 4% of the rated flow as shown in Figure 14. The total power continues to decrease due to the negative net reactivity. The reactor enters a long-term cooling phase.Figure 21
Row 2 and row 6 PIOTAs outlet temperature.Figure 22
Row 2 and row 6 PIOTAs outlet temperature (0 s–200 s).Figure 23
Dimensionless power-flow ratio.In general, the predicted results of SAC-3D are in good agreement with the trend of the experimental data, but there are still deviations in the temperature values. The outlet temperature of the PIOTA in row 2 is consistently lower than the measured value, which is caused by the bias of the simulated value of the total power transient variation. Similarly, the simulated outlet temperature of the PIOTA in row 6 is lower than the measured data until about 20 s but starts to be higher than the measured data after then. There is a maximum deviation of about 40 degrees Celsius at the peak. Unlike the PIOTA in row 2, the PIOTA in row 6 is located at the outermost ring of the fuel assembly and is next to a GEM assembly and a reflector assembly. We believe that there are at least two possible reasons for this anomaly. One is that even though the transient's overall mass flow rate is in good agreement with the measurement data, it is still possible that the simulated mass flow rate variation within each channel may not match the actual situation all that well. The other important factor is the fact that the channel-to-channel heat transfer is not considered in the present work. The neglect of heat transfer between the PIOTA in row 6 and GEM and reflector assembly introduced positive temperature errors.


## 5. Conclusions and Perspective

In the blind calculation phase of the IAEA FFTF LOFWOS test #13 CRP, we modeled the FFTF reactor system by using the SAC-3D system analysis code. The SAC-3D calculation results have reasonably good agreement with the measured data. After benchmark analyzing the EBR-II [2] and FFTF reactors, SAC-3D shows its good capability to predict the fast reactor behaviors under ATWS transient conditions. However, the SAC-3D and FFTF models still have room for improvement. The deviations in total power transient calculations suggest further improvement in the reactivity feedback models as well as the decay heat prediction model. Anomalies in the transient temperature prediction of the PIOTA in row 6 indicate that the thermal-hydraulics calculation module needs to consider the heat transfer between channels. In the next phase of the work, we will examine the deviations between simulation results further and refine the FFTF model based on the blind phase benchmark analysis results. The specific calculation modules of SAC-3D will also be improved to predict the fast reactor's behavior better.



---

*Source: 5843910-2021-11-10.xml*