Effect of heat transfer correlations on the fuel temperature prediction of SCWRs

In this paper, we present a numerical analysis of the effect of different heat transfer correlations on the prediction of the cladding wall temperature in a supercritical water reactor at nominal operating conditions. The neutronics process with temperature feedback effects, the heat transfer in the fuel rod, and the thermalhydraulics in the core were simulated with a three-pass core design.


Introduction
The super critical water reactor (SCWR) is one of the most promising and innovative designs selected by the Generation IV International Forum.This is a very high-pressure water-cooled reactor which will operate at conditions above the thermodynamic critical point.Water enters the reactor core and then exits without change of phase, i.e., no water/steam separation is necessary.There is an increase of thermal efficiency of current nuclear power plants from 30-35% to approximately 45-50%.
Figure 1 shows the difference in the operating conditions of current generation reactor systems in comparison to SCWRs.Compared to existing pressurized water reactors (PWRs), in SCWRs the target is to increase the coolant pressure from 10-16 MPa to about 25 MPa; the inlet temperature to about 350 °C, and the outlet temperature to about 625 °C [1].
In this paper, we presented a numerical analysis of the effect of different heat transfer correlations on the prediction of fuel and wall cladding temperatures in a supercritical water reactor.The neutronics process with temperature feedback effects, the heat transfer in the fuel rod and the thermal-hydraulics in the core were simulated.Special attention was given to the thermal-hydraulics, which uses a three-pass core design with multiple heat-up steps, where each step was simulated using an average channel.The first pass called "evaporator" is located in the center of the core.In this region, the moderator water flows downward in gaps between assembly boxes and inside the moderator tubes.The moderator water, heated-up through its path downward to the lower plenum, is mixed with the coolant coming from the downcomer reaching an inlet temperature of around 583 K.The evaporator heats the coolant up to 663 K, flowing upward and around the fuel rods, resulting in an outlet temperature 5 K higher than the pseudo-critical temperature of 557.7 K at a pressure of 25 MPa.The second pass, called "superheater", with downward flow, heats the coolant up to 706 K.After a second mixing in an outer mixing plenum below the core, the coolant will finally be heated up to 803 K with an upward flow in a second superheater (the third pass) located at the core periphery.A transient one-dimensional radial conduction model was applied in the fuel rod for each cell in the axial coordinate.Energy balances for the coolant have been implemented using a steady state and a onedimensional model for the axial coordinate.Fuel lattice neutronics calculations were performed with the HELIOS-2 code and the reactivity coefficients were used to evaluate the reactivity effects due to changes in the fuel temperature and in the supercritical water density for 177 energy groups.Due to the strong variation of coolant density through the core, five densities were considered.This safety parameter is calculated in order to evaluate the variation of the reactivity due to the Doppler effect, as a function of the fuel temperature, which is related to the resonances broadening when the fuel temperature increases.The coupling of neutronics with the heat transfer in the fuel rod, and the thermal-hydraulics is presented, and numerical experiments due to changes in the mass flow rate were accomplished in this study.Effects on fuel temperature predictions with improved heat transfer correlations and classical heat transfer correlations were also compared.

Supercritical fluids
The behavior of liquid and gas density with pressure and temperature is illustrated in Figure 2. When the pressure and temperatures are low, there is a significant density difference between the liquid and the gas states.Near the critical point, the density difference between the liquid and gas is small, and above the critical point, the densities of the liquid and the gas have become equal.The heat transfer process, at critical and supercritical pressures, is influenced by the significant changes in thermophysical properties, as is observed in Figure 3 for specific heat, thermal conductivity, and density obtained from thermal properties taken from [3].The most significant thermophysical property variations occur near the critical and pseudocritical points.For example, the specific heat of water has a maximum value at the critical point.The exact temperature that corresponds to the specific heat peak at pressures above the critical pressure is known as the pseudocritical temperature [4].

Supercritical water heat-transfer correlation
The practical prediction methods for heat transfer at supercritical pressures are presented in [1,4].The supercritical water heat transfer correlations applied in this work are shown in Table 1.Dimensionless numbers used in Table 1 are given by: which is the Nussel number.Here H ∞ is the heat transfer coefficient, k is the thermal conductivity and D H is the hydraulic characteristic length.The subscript b means that bulk-fluid temperature is used to calculate the thermophysical properties.These properties can also be calculated with the wall temperature, which will be specified with a subscript w.The Reynolds number is defined by: where G is the mass flux and m is the viscosity.The Prandtl number is defined as: where Cp is the specific heat.The heat transfer coefficient is used in the boundary condition given below in equation (6), and H ∞ represents the heat transfer from the wall to the coolant.McAdams [6] proposed the use of the Dittus-Boelter correlation for forced-convective heat transfer in turbulent flows at subcritical pressures.The only difference between the Dittus-Boelter and McAdams correlations is that the latter has a larger coefficient.According to Schnurr et al. [10], it agrees with experimental data.However, it was noted that the correlation might produce unrealistic temperature results near the critical and pseudocritical points, due to it being very sensitive to variations in the thermophysical properties.Bishop et al. [7] conducted experiments in supercritical water flowing upward inside bare tubes and annuli, within the following range of operating parameters: P = 22.8-27.6MPa, T b = 282-527 °C, m = 651--3662 kg/m 2 s and q = 0.31-3.46MW/m 2 .Their data for heat transfer in tubes were generalized with a fit of ±15%.This correlation uses a cross-sectional averaged Prandtl number and the final term in the correlation (1 + 2.4 D/x) accounts for the entrance-region effect.Bishop et al.'s correlation was modified and used without the entrance-region term, because this term depends significantly on the particular design of the inlet of the bare test section.
Swenson et al. [8] have suggested a correlation in which thermophysical properties are mainly based on a wall temperature, as they found that conventional correlations, which use a bulk-fluid temperature as a basis for calculating the majority of thermophysical properties, did not work as well.
A dimensional analysis was performed by Mokry et al. [9] in order to obtain a general empirical form of correlation for the heat transfer calculations, and as a result of the experimental data analysis, two correlations for the heat transfer coefficient at supercritical water conditions were obtained.
In the core layout of the SCWR under study, water, as the working fluid, is guided three times through the core (twice up and once down).This design is called the threepass core concept.The first pass, called the evaporator, is situated in the center of the core.In this region, the moderator water flows downward in gaps between assembly boxes and inside the moderator tubes.The moderator water, in its downward path to the lower plenum is heated up, and is mixed with the coolant (1200 kg/s as inlet mass flow) which comes from the downcomer, thereby reaching an inlet temperature of around 583 K.The evaporator heats the coolant up to 663 K, flowing upward around the fuel rods, resulting in an outlet temperature 5 K higher than the pseudo-critical temperature of 557.7 K at a pressure of 25 MPa.An inner steam plenum above the core eliminates hot streaks.The second pass, called superheater, with a downward flow, heats the coolant up to 706 K.After a second mixing in an outer mixing plenum below the core, the coolant is finally heated up to 803K with an upward flow in a second superheater located at the core periphery, known as the third pass.Each pass, the evaporator and both superheaters, is built of 52 fuel assembly clusters as shown in Figure 4 [11].Therefore the complete reactor core is composed of 156 assembly clusters.
The fuel assembly design is taken from the European high performance light water reactor (HPLWR) concept.The 7 Â 7 square arrangement design, with 40 fuel rods distributed in dual rows, and a single water tube replacing 9 fuel rods was used [12].The fuel rods and the water tube are housed within the assembly box and grouped in a cluster of 9 assemblies, in a 3 Â 3 arrangement with similar dimensions to a PWR assembly.As found in the PWR, control rods are inserted from the core top into 5 of the 9 water tubes of a cluster (Fig. 4b).The structural material for cladding, assembly boxes and water tubes is stainless steel.The main reactor parameters are presented in Table 2.

Implementation of the heat transfer correlations
In order to analyze the effect of different heat transfer correlations on the prediction of the wall temperature of the fuel rods, the SCWR numerical code developed by Barragán-Martínez [15] was applied using the HTCs shown in Table 1.The numerical model of the heat transfer processes in the fuel element of the HPLWR was obtained using the numerical model of typical reactors [16].The supercritical water reactor is integrated of cylindrical fuel elements which contain ceramic pellets inside the cladding.
Then, the effect of heat transfer correlations on the fuel temperature prediction of SCWRs was conducted with numerical experiments.

Fuel heat transfer model
A detailed multi-node fuel pin model was developed for this study.The fuel heat transfer formulation is based on the following fundamental assumptions: axis-symmetric radial heat transfer, the heat conduction in the axial direction is negligible, the volumetric heat rate generation in the fuel is uniform in each radial node, and storage of heat in the fuel cladding and gap is negligible.Under these assumptions, the transient temperature distribution in the fuel pin, and the initial and boundary conditions are given in the following conditions: ∂T ∂r þ q 000 ðtÞ; at r r r f ; ð4Þ In equation (4) q 000 (t) = 0, for r f r cl .In these equations, r is the cylindrical radial coordinate, r 0 , r f and r cl are the centroid, fuel and clad radius, respectively, q (t) = P(t)/V f at each axial node, where P is the neutronics power, T m is the moderator temperature, and H ∞ is the convective heat transfer coefficient.
The differential equations described previously are transformed into discrete equations using the control volume formulation technique in an implicit form [17].The control volume formulation enables the equations for fuel, gap, and cladding to be written as a single set of algebraic equations for the sweep in the radial direction:  where T tþDt jÀ1 , T tþDt j and T tþDt jþ1 are unknowns, a j , b j , c j and d j are coefficients, which are computed at the time t.When these equations are put into a matrix form, the coefficient matrix is tridiagonal.The solution procedure for the tridiagonal system is the Thomas algorithm, which is the most efficient algorithm for this type of matrices.The coefficients a j , b j , and c j are dependent on thermophysical properties, i.e., thermal conductivity, density and specific heat; and since they are function of T tþDt j , at least one iteration is needed.

Thermal-hydraulics model
The basic equations for describing the thermal-hydraulics behavior in the three representative heated channels (one channel for each pass core) assuming the supercritical fluid is a single phase fluid, are presented as following.Incompressible flow was also considered in this study, i.e., the mass flux (G) is a constant.Under this consideration, the energy equation at steady state is shown as follows: where T is the temperature, f is the friction factor, P H is the heated perimeter, A f is the flow area.The heat transfer from the wall to the coolant is obtained with Newton's law of cooling: The temperature in each node of the channel is obtained numerically as: where Dz is the node length and i is the node number.

Reactor power model
The reactor power is given by where F(z) is the axial power factor, P 0 is nominal power and n(t) is the normalized neutron flux, which is calculated by using a point reactor kinetics model with six groups of delayed neutrons: where C i is a delayed neutron concentration of the ith precursor group normalized with the steady-state neutron density, r is the net reactivity, b is the neutron delay fraction, L is the neutron generation time and b i is the portion of neutrons generated by the ith group.The initial conditions are given by n(0) = n 0 and c i (0) = b i n 0 /Ll i .The parameters of the kinetics model are presented in Table 3.
The net reactivity in this work includes three main components: Doppler effects due to fuel temperature, coolant density, and reactor control rods.
The kinetics point equations are stiff in the coefficients because they differ in several orders of magnitude.The implicit variable integration method was used to solve equation ( 13), and the Euler method in an explicit form was used to solve the delayed precursor concentration given by equation (14).
The reactivity coefficient due to variations in fuel temperature was studied for the square fuel assembly design proposed by [18].The calculations were done for the fuel assembly model along the active core height.Due to the strong variation of coolant density through the axial direction of the core, five densities: 0.74, 0.45, 0.31, 0.17 and 0.09 g/cm 3 were considered.This safety parameter is calculated in order to evaluate the variation of the reactivity due to the Doppler Effect, as a function of the fuel temperature, which is related to the resonances broadening when the temperature increases.The values of the reactivity as a function of the coolant density and fuel temperature are presented in Figure 5.The values of the infinite multiplication factor obtained with HELIOS-2 for 177 energy groups were used to determine the reactivity [19].

Representative SCWR nodalization
The fuel rod temperature distribution was obtained for the radial nodes at each of the twenty one thermal-hydraulics axial nodes in the core.The arrangement of the computational nodes of the thermal-hydraulics model is illustrated in Figure 6.
Figure 7 shows the grid used in calculations.Half control volume near the boundary, radial nodes 1, 2, 3, 4, and 5 for the fuel; radial node 6 was used for the gap; radial nodes 7 and 8 for the clad.Radial nodes 1 and 8 were used for the boundary condition.

Numerical experiments
Each channel in the core was based on an hydraulic unit cell whose parameters are: P H = 0.025 m, D H = 0.054 m, and A f = 0.34 m 2 .The parameters of the fuel element are:     1.It should be noted that the last node temperature (at 4 m) is practically the same, and the trend is very similar for all the correlations, except for a short zone where the Swenson correlation yields a lower temperature while Mokry's correlations (both preliminary and final) yield a higher temperature, the same was noted for the Bishop's correlations (with and without ERE).
In Figure 9 the results for Channel 2 are presented, showing the wall temperature behavior for the correlations in Table 1.Similar results were obtained, however contrary to what was observed in Channel 1, the Swenson's correlation yields slightly higher temperatures along the entire channel meanwhile the Bishop's (with and without ERE) and Mokry's correlations yield slightly lower temperatures along the entire channel.
Figure 10 presents the results for Channel 3, showing the wall temperature behavior for the correlations presented in Table 1.In this case, the trends that most resemble each other are presented.Again, the Swenson's correlation deviates the most, yielding slightly higher temperatures than other correlations.
In Figure 11 the results along the three channels are presented.It should be noted that Swenson's correlation is the one with greater deviation from Dittus-Boelter's correlation, with a difference of 10 K in Channel 1.There is a wall temperature reduction at the end of each channel; especially for Channel 2, and this is due to the axial distribution of thermal power which has a minimum in this bottom core zone.This is an undesired result of the three pass core concept.Another finding in this numerical analysis was that Swenson's correlation gave the most conservative predictions, in terms of safety, because higher temperatures are calculated due to the use of the wall temperature for the Re and Pr calculations, while the other correlations use the bulk temperature.

Conclusions
The correlation, which agrees most with Dittus-Boelter, is McAdams.The only difference in the equation is the value of the coefficient.Bishop's correlations, with and without Entrance-Region Effect (ERE) have little differences among them in the prediction of the wall temperatures, meaning that, for this simulation the ERE is not important; predictions compared to the Dittus-Boelter correlation are a little higher in the first channel and slightly lower in Channels 2 and 3.With preliminary and final Mokry's correlations, higher temperature predictions were found in Channel 1, but were very similar to Dittus-Boelter in Channels 2 and 3. Swenson's correlation showed the most deviated results, yielding lower temperatures in the first channel and higher in Channels 2 and 3.
Swenson's correlation uses the wall temperature for calculating the Re and Pr numbers, while the others used the bulk temperature and we found the greatest differences compared to other HTCs.For this reason Swenson's correlation could be very useful in order to find the most conservative results for Channel 3, where high wall temperatures could affect the fuel rod integrity.

Fig. 3 .
Fig. 3. Behavior of the specific heat (C p ), thermal conductivity (k) and density (r), as a function of the temperature at 25 MPa.

Fig. 6 .
Fig.6.Arrangement of the computational nodes in the thermal-hydraulics core model of the SCWR.

Figure 8
presents the results for Channel 1, showing the wall temperature behavior for different correlations presented in Table

Fig. 7 .
Fig. 7. Arrangement of the computational cells of fuel, gap, and clad.

Fig. 10 .
Fig. 10.Simulation results for Channel 3 showing the Wall Temperature behavior for different HTCs.
Special thanks to the National Council for Sciences and Technology (CONACYT) for the scholarship provided to the Master Student Erick G. Espinosa-Martinez, and to the National Autonomous University of Mexico for the PAPIIT IN113213 project funds.