Issue
EPJ Nuclear Sci. Technol.
Volume 5, 2019
Progress in the Science and Technology of Nuclear Reactors using Molten Salts
Article Number 16
Number of page(s) 13
Section Physics
DOI https://doi.org/10.1051/epjn/2019027
Published online 14 November 2019

© R. Freile and M. Kimber, published by EDP Sciences, 2019

Licence Creative Commons
This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1 Introduction and background

During the last 50 years there has been a rising interest on the use of molten salts as a heat transfer fluid. In particular, they have been considered among the best candidates for the advanced reactor designs of the Generation IV reactor concepts [1].

Molten salts’ primary function in a nuclear reactor is to act as a coolant, extracting the heat resulting from the nuclear fission. These designs may have variations according to how the fuel is arranged in the core. For example, the fuel may be present in the form of a ceramic fuel in a prism or pebble bed core design (Liquid-salt-very-high-temperature reactor) or rather dissolved in the salt itself (i.e. liquid fuel). One very promising design is the liquid-fueled molten salt reactor, which is the only design of the Generation IV concepts which employs a liquid fuel. Fluorides of fissile and/or fertile elements such as UF4, PuF3 and/or ThF4 are combined with carrier salts to form fluids. The most famous carrier salt candidates are FLiNaK (46.5 LiF-11.5 NaF-42 KF mol%) and FLiBe (66 LiF-34 BeF2 mol%).

In a liquid fuel molten salt reactor (MSR) a key factor towards its modeling is the strong coupling between the physics describing the neutronics and the thermal-hydraulics. The velocity field affects the neutron precursors transport and the temperature field affects the neutron population through Doppler effects and density changes in the coolant. Hence, it is of paramount importance to be able to describe the fluid flow and the thermal characteristics of the molten salt fluid. In order to do so, it is required that the turbulence models used are valid for the molten salt flow and that the FLiNaK material properties are accurately described.

Historically, the investigation of molten salts began in the late 1950s at Oak Ridge National Laboratory, aimed at the development of a molten-salt nuclear reactor to power an airplane. In 1954 the first molten salt reactor for the Aircraft Reactor Experiment was built at Oak Ridge to investigate the use of molten fluoride fuels for an aircraft propulsion reactor.

During the same time period, forced convection experiments were performed to examine heat transfer properties of molten salts, in particular FLiNaK. In 1954 Grele and Gedeon [2] and Hoffman and Lones [3] performed experimental measurements using FLiNaK flowing through a heated test section with different piping materials. Nusselt numbers were calculated and contrasted with the Dittus-Boelter correlation, depicting an underprediction of about 50% in value. Final results revealed good agreement when using Ni and SS316, but as much as 50% of underprediction using Inconel. The results obtained by these experimentalists is shown in Figure 1.

Later inthe 1970’s, Vriesema performed measurements on forced convection in a vertical pipe [4]. The results of Nusselt number were 15% off the Dittus Boelter prediction, as shown in Figure 2.

Ignat’ ev et al. [5] also performed heat transfer studies using FLiNaK in a circular tube of iron-based steel. The results were compared with the Petukhov-Kirillov formula, yielding a maximum error of 10% with respect to their experimental data.

Throughout all these experiments appearing in the literature, they all had the common problem of not knowing precisely the value of the thermophysical properties of the molten salt. In particular, the wide range of values adopted by these researchers for thermal conductivity is worth noting. Grele and Gedeon and Hoffman and Loans used for their Nusselt calculations a thermal conductivity of λ = 4.5 W/mK, while Vriesema used λ = 1.3 W/mK. In the Greleand Gedeon report it was mentioned that the reason for the huge discrepancy between their measurements and Dittus Boelter was due to the formation of a film associated with the chromium that was corroded from the Inconel piping, resulting in an additional resistance to heat transfer in the Inconel tube. Studies done by [8] suggest thatanother possible reason for the existent large discrepancy was the adopted value of Flinak thermal conductivity.

Experimentally determining thermal conductivity for molten salts is a challenging task, as it is necessary to differentiate the contributions from the radiation, natural convection and conduction heat transfer modes in a given experiment [6, 8]. Because the thermal conductivity values in the past overestimated by a factor of 4 the most widely accepted values nowadays [7], it can be argued that all the heat transfer modes were not taken into account during the thermal conductivity measurements in that era.

In 1962, Ewing [6] was the first to realize that the thermal radiation in molten salts could significantly affect the experimental measurements for the conductivity. He described the measured value as an effective thermal conductivity, accounting for a molecular conduction plus a thermal radiation conduction. Ewing concluded that the effective thermal conductivity values were very sensitive to the conditions of the experiment. In his measurements, the FLiNaK conductivity values ranged from 0.6 W/m K to 5 W/m K. The former was constant during each run equal to λ = 0.6 W/m K.

Another set of experiments found in literature for measuring the thermal conductivity of FLiNaK was done by Smirnov [7] in 1987. He used the method of coaxial cylinders made of platinum, with radiation heat transfer taken into account. He proposed an empirical correlation for the thermal conductivity equal to λ = 0.36 + 0.00056.T(W/m K) over the range of 517 °C to 817°C.

In 2009, Ambrosek [8] reanalyzed each of the previously mentioned experiments using the most widely accepted thermophysical properties of FLiNaK. In particular, he used Smirnov’s correlation, and found thermal conductivity values of k = 0.81 − 0.93 W/m K over the temperature range of validity. The final result of this analysis agreed within 15% with respect to Dittus Boelter correlation for Inconel piping. The reanalyzed data by Ambrosek is shown in Figure 3. Even though in this figure, the results for Ni and SS316 pipes was not shown, it must be said that they presented a large discrepancy when compared with the Dittus Boelter correlation. No explanation could be found for these results.

After this work done by Ambrosek, it was concluded that the main reason why the first experimentalists had observed results in great discordance with the known correlations was the misprediction in the thermal conductivity of FLiNaK used for their calculations.

In 2013, Idaho National Laboratory published an Engineering database of liquid salt thermophysical properties [9]. It consists of acompilation of the different commonly accepted correlations for the temperature dependent FLiNaK properties such as specific heat, density, viscosity and thermal conductivity. In Table 1 some of the most widely accepted correlations for these properties are shown. Note that for the thermal conductivity there is only one expression in the table, which is the most widely accepted temperature dependent expression. However, according to studies done by Williams et al. [12], the thermal conductivity value should range between 0.6 W/m K and 1.0 W/m K at 973 K.

As noted earlier, it is unclear which set of expressions would be best to use. With this in mind, the authors in [9] determined the effect of uncertainty in the resulting Nusselt number, but assumed complete validity of the Dittus Boelter equation. For the worst case scenario, the maximum reported error due to thermophysical uncertainties in the Nusselt number was + ∕−8%. However, those authors suggested a more detailed analysis should be conducted, one that does not rely on Nusselt number correlations or the assumption of constant properties.

In the scientific community other efforts were made to perform sensitivity analysis of heat transfer in molten salts. Xu et al. [10] conducted a sensitivity study of heat transfer mechanisms in a packed-bed molten salt thermocline thermal storage system using a two-phase model. Different interstitial heat transfer correlations, which are functions of the Reynolds number and the porosity, were used to compare the thermal performance of the storage system. Sabharwall et al. [11] performed hand-made calculations of the sensitivity coefficients for the different thermophysical properties in the VHTR (very high temperature reactor). In their calculations, an analytical equation for the heat transfer coefficient was obtained by using the Dittus Boelter equation and using a fixed pressure drop with friction coefficients for the velocity estimation, which was made using the Blasius formula. Partial derivatives of the heat transfer coefficient with respect to each individual thermophysical property were taken and dimensionless sensitivity coefficients were obtained. Therefore, in the calculation of the sensitivity coefficients, they assumed that Dittus Boelter is valid and that the bulk temperature of the salt does not strongly depend on the molten salt properties.

To the best of our knowledge, there is no definitive expression that best defines molten salt material properties, and more importantly,a thorough investigation is needed to truly quantify the impact uncertainty in material properties has on the pertinent thermal hydraulics. Several expressions for each material property are found in the literature, with noticeable differences in the values one would predict. In the present work, computational fluid dynamics (CFD) was used to obtain sensitivity coefficients for each material property. With this information, the effect of the uncertainties in the material properties on the heat transfer coefficient was accounted for. The CFD simulations using RANS models were compared to a published experiment on molten salt done by Grele and Gedeon [2] and the adequacy of different Nusselt number correlations was assessed.

thumbnail Fig. 1

Experimental results obtained by Grele and Gedeon [2] and Hoffman and Loans [3] compared with Dittus Boelter correlation.

thumbnail Fig. 2

Experimental results obtained by Vriesema [4] compared with Dittus Boelter.

thumbnail Fig. 3

Grele and Gedeon, Hoffman and Loans and Vriesema results for Inconel piping reanalyzed by Ambrosek et al. [8].

Table 1

Correlations for thermophysical properties of FLiNaK taken from [9].

2 Model

2.1 Hydrodynamics

In this work, equations of flow motion were solved using the Reynolds Averaged Navier Stokes models. Using this approach, the instantaneous velocity is expressed as the sum of a mean velocity , and a fluctuating part, , such that: (1)

where

Replacing the instantaneous velocity in the Navier Stokes equations, and time averaging both yields the Reynolds averaged equations ofmotion: (2) (3)

It must be clarified that the compressibility effect (change in the density of the fluid) can be caused either by temperature or by pressure. In the case of molten salt flow in a pipe, from the standpoint of the pressure influence on the density, it can be said that it will remain constant. The way to justify this statement is by analyzing the Mach number, which in this particular case is much lower than one. On the other hand, in heated flows, there will be a density variation due to its dependence on temperature. Although the dilatation term was modeled, it can be said that its value will still be negligible compared with the rate of strain tensor.

The last term on the right hand side of equation (3) represents the effects of turbulence and it is commonly known as the Reynolds stress tensor. Note that for this Reynolds stresses an overbar is used to denote the time averaged quantities of those variables. In this work, the Boussinesq approximation was used. The expression for the Reynolds tensor using this approximation is shown in the next equation: (4)

where k is defined as the turbulent kinetic energy and μt is the turbulent viscosity, whose expression varies accordingly to the model used.

In this work a two-equation model was used. In these types of models the turbulent viscosity μt is taken as a function of the turbulent kinetic energy and the turbulent dissipation . Particularly for this work, the k-ω SST model was used, where is called specific dissipation or turbulent frequency.

The shear-stress transport (SST) k-ω model was developed by Menter [15] with the objective of combining two of the most popular two-equation models, k-ω and k-ɛ, through the use of blending functions. Both of the aforementioned models are multiplied by a blending function and both models are added together. The blending function is generally a hyperbolic tangent whose range is restricted between zero and one. It is designed in such a way that near the wall, the standard k-ω model is predominant, and zero away from the surface, which activates the k-ɛ model in the free stream region. By doing this, the k-ω SST model is directly usable all the way down to the wall through the viscous sub-layer, hence being able to resolve turbulent parameters up to the wall region. In addition, in the free stream region it avoids the common k-ω problem of being overly sensitive to the inlet free stream turbulence properties and takes advantage of the k-ɛ model’s advantage for free-shear flows.

There are several works which have been done on the CFD RANS modeling of molten salts. Ferng et al. [20] proposed a CFD methodology for investigating thermal-hydraulic characteristics of FLiNaK in a pipe geometry using kɛ model. Results were compared to existing correlations. Chen et al. [19] performed calculations of Nusselt numbers for Hitec molten salt using four different RANS models and compared them with a present experiment of Hitec and several other experiments with different salts.

This turbulence model was chosen because of the good performance it has shown in other simulations related to the geometry of interest (i.e., pipe flow or channel flow). Kim et al. [16] investigated the effects of non-uniformity of fluid properties on forced convection heat transfer using different turbulence models by comparing against experiments. In this study, the effect of the non-uniformity of fluid properties was accounted for by applying a factor to the value of the Nusselt number for constant properties. The prediction with the kω model formulation was in good agreement with the experimental results. Menter et al. [17] compared the performance of the kω SST model in heat transfer applications with other turbulence models such as Low-Re models and k-epsilon models against experiments for different geometries. The best overall performance was achieved with the kω SST model. Given the conclusions of these two studies, the kω SST model is employed in the current work as well.

2.2 Heat transfer

Similarly to what was shown in equation 1, the energy E and enthalpy H are also decomposed into a mean and a fluctuating value. In order to compute the temperature distribution in the fluid (Eq. (5)) and in the solid (Eq. (6)), the equations of energy conservation were used: (5) (6)

where , and is defined as a volumetric heat source in the solid. The term is commonly known as the turbulent heat flux.

In order to model this last term, a similar concept to the one used in Boussinesq approximation (Eq. (4)) is used: (7)

where αt is the eddy-diffusivity and Prt is the turbulent Prandtl number. This number was given a generally used constant value of 0.85 for the simulations [15].

Regarding radiation effects, results provided by Ambrosek et al. [8] proved that, because of the transparency of FLiNaK [6], the amount of energy transferred by radiation can be significant in applications involving high temperatures (T = 1123 K) and laminar flow conditions (Re < 500) in pipes with a diameter of 1 cm or greater. In the current experimental configuration, the pipe diameter is too small to take into account radiation heat transfer, as the convective heat transfer coefficient prevails over the radiative one. This was also proven by Ethan S. Chaleff et al. [22].

3 Methodology

3.1 Computational domain

RANS simulations were carried out in this study using ANSYS Fluent v19 commercial package. For this work the k-ω SST model was used for reasons previously noted. CFD codes require generating a computational domain and mesh to simulate flow behaviors in a specific geometry. The computational domain modeling the experimental test section in Grele and Gedeon experiment [2] was completed using ICEM CFD 19.0. The computational domain is shown in Figure 4.

The model consisted of a solid volume modeling the Inconel X pipe wall and a fluid volume to represent the molten salt flow. The diameters, the pipe wall thickness and the test section’s length dimensions were in agreement with the test section used in [2]. Specifically, it had an outside diameter of 3/8 inches and a wall thickness of 0.065 inches. The total length was 24 inches.

A mesh convergence study using Richardson extrapolation [18] was made to assure that the model was producing a mathematically accurate solution. The number of elements was gradually increased until the wall temperature error for the highest Reynolds number casewas below 1E-2 K. The final total number of elements in the mesh was 1,054,585. The minimum orthogonal quality was 5.14e−1 and the maximum aspect ratio was 50.4 in the cells which are the nearest to the wall. As seen in Figure 4, in the radial direction a 30 × 30 uniform square cartesian grid was used at the center of the pipe. The vertices of the square were equidistant from the origin at a distance of 1.8 mm. Thirty additional radial layers of cells exist between the bounds of the square and the inner diameter of the pipe. A growth ratio of 1.005 from the wall was imposed in this section, and the first mesh division was set at 2E-6 m. This refinement in the fluid region near the wall was made so that the model could resolve the equations for the temperatures and velocities up to the viscous sub-layer. In internal flows, the distance from the wall is generally represented by using the y+ coordinate, which is defined as , where y is the distance from the wall and τw is the wall shear stress. For the highest Reynolds case, the mesh contained at least one grid point below normalized wall distance, i.e y+ =1, and at least five grid points in the viscous sublayer region, i.e. y+ = 5. In the solid region, 8 equidistant divisions were imposed radially. In the axial direction, two hundred uniformly spaced divisions were imposed.

The numerical solution was obtained using the SIMPLE algorithm [21]. The gradients were computed by using the Green-Gauss Cell-Based method, which is taken from an average of the values at the neighboring cell centers in a certain direction. For the convection terms, face values of different variables are needed. For the momentum, turbulent kinetic energy, dissipation and energy equation a second order upwind scheme was used, which takes information from the cell-center value and gradient value in the upstream cell. The simulations were considered converged once the residuals levels became lower than 1E-8.

thumbnail Fig. 4

Perpendicular section of the computational domain used in the simulations.

3.2 Boundary conditions

The boundary conditions applied to the numerical model include standard hydrodynamic conditions for pipe flow as well as thermal boundary conditions in both the fluid and solid regions of the model. The hydrodynamic boundary conditions are provided in Table 2, and include a prescribed hydrodynamically-developed velocity at the inlet, pressure outlet and no slip conditions at the wall. The inlet velocity is reported in [2], and was experimentally determined for each trial by measuring the amount of time it took for the fluid level to rise from one level to another in a volume-measuring tank locate d within the loop. Those authors used a pair of contact points and an electric stopwatch to make this calculation. The turbulence inlet boundary conditions were inputted by using standard expressions for a fully developed internal flow, defining turbulence intensity and a reference length.

Solving the energy equation also required appropriate boundary conditions to be defined and are needed for both the fluid and the solid within the domain. These are shown in Tables 3 and 4, respectively. The fluid inlet temperature for each trial is taken from [2] and the outlet is considered adiabatic. At the surface between the molten salt and the solid tube wall, an interface condition is applied which serves to thermally couple the two element types (i.e., solid and fluid). For the solid, this interface condition is also applied, while the remaining solid surfaces are all adiabatic. The coupled boundary condition was set in the interface. The solution for the temperatures in the interface are solved for by using the information of the temperatures of the interface adjacent cells, iterating until the energy balance is satisfied at the interface.

Enhanced wall treatment was used during the simulations, which means that for the dimensionless velocity, a blending function is used to smoothly merge both the log layer and the viscous layer. For the thermal formulation, an elliptic blending is also used for merging the laminar and logarithmic profiles accordingly.

In Grele and Gedeon experiment, electrical heaters were fixed around the pipe circumference which provided the required heating. To model this, cell zone conditions were imposed in the solid region in order to add the volumetric heat source due to electrical heating. The constant volumetric heat source may have discrepancies with the actual experiment because of several factors such as non-uniformity of the electrical current and because of the variation of the pipe temperature. However, because of the lack of data regarding the heat flux profiles, it was decided that the constant volumetric heat source was the best choice for reproducing the experiment as it represents the simplest solution for applying the thermal load.

In each one of the cases, the total heat supplied to the fluid was computed using the temperatures at the inlet and the outlet together with the choice of material properties (Tab. 1), particularly density and specific heat. The formula used was the following: (8)

where superscripts x and y refer to the different possibilities of thermophysical properties (Tab. 1), Tb is the bulk temperature calculated as and A is the cross- sectional area of the pipe. After the total heat supplied to the fluid was computed, using the assumption of uniformly distributed heat, the volumetric heat source was inputted as a cell zone condition in the solid region.

Table 2

Hydrodynamic boundary conditions used.

Table 3

Energy boundary conditions used for the fluid region.

Table 4

Energy boundary conditions used for the solid region.

4 Results

In order to assess the RANS models for molten salt flows (FLiNaK), different cases in the Grele and Gedeon experiment were simulated and compared to their experimental results. In the experiment, the external average wall temperature was reported in each run, which was found by simply averaging the 12 thermocouples located on the pipe wall at different axial positions. The wall temperature in the external wall of the modeled solid pipe at each one of the thermocouples’ positions was averaged in order to directly compare with the reported external wall temperature from the experiment.

The thermophysical properties of FLiNaK were considered as temperature dependent in the simulations, implemented as user defined functions in Fluent. Due to the lack of trustworthy material properties, for each case, different combinations of the material properties of FLiNaK presented in Table 1 were used. The thermal conductivity used was given by Smirnov’s correlation, as the majority of the published literature consider this correlation as the most accurate for FLiNaK. In regards to the solid region, the Inconel X pipe thermophysical properties were taken from [23].

From the 52 cases presented by Grele and Gedeon, 10 specific cases were simulated, and encompass the range of operating conditions encountered by the experimentalists. The key parameters taken from the experimental data were:

  • for the fluid: the velocity inlet, the temperature inlet;

  • for the solid: using equation (8), taking the inlet and outlet temperature, the uniform volumetric heat flux was calculated and imposed as a cell zone condition.

The 10 cases selected for the simulations are shown in Table 5.

For clarification, the range for the Reynolds and Prandtl numbers observed in column 6 and 7 of Table 5 corresponded to the minimum and maximum values using the aforementioned different combination of thermophysical properties.

Taking into account the two correlations for density, the two correlations for specific heat and the four correlations for viscosity, the total amount of 16 simulations were run per case, yielding a total of 160 simulations. Thus, in order to proceed in a practical way, a loop which went through all of the possible combinations was programmed in Scheme language, a language seamlessly interpreted by Fluent. In this script, the case was initialized, the properties were changed accordingly, the simulation was run andfinally results were exported. Using Matlab, a table containing all the pertinent results was created automatically. In particular, this code calculated the heat transfer coefficient and compared it with Dittus-Boelter and Sieder Tate correlation, both of which useReynolds and Prandtl numbers as inputs, and therefore a spread exists in that data as well due to dependence on thermophysical properties. Unsurprisingly, the Nusselt number spread using the Sieder Tate correlation was smaller in comparison to that using the Dittus Boelter correlation. On average, the simulation results showed a 1.4% difference when compared to the Sieder Tate correlation and a 4.9% of the Dittus Boelter correlation. This showed that the variation of the viscosity values at the wall surface with respect to the bulk viscosity was significant and could not be neglected.

The results for the averaged outer wall temperature are shown in Figure 5. On the y-axis, the difference between the experimental results and the simulation results for the averaged wall temperature was plotted. The x-axis denotes the specific case number as previously noted.

The most noticeable conclusion from Figure 5 is that the variations on the predicted thermophysical properties can lead to a wall temperature spread of almost 15 ℃. The constant line with a zero value represents the match between the simulation and the experimental results.

The initial thought with running these experiments was to identify the set of material properties that best describe the experimental results. The properties’ correlations that best describe the experiment were not consistent between any two experiments, meaning that the best set of material properties could not be found. It should be noted that even though for some cases (e.g., 30 and 46), results suggest that there is no overlap between experimental results and simulation results, the reader should bear in mind that experimental uncertainties were not reported and cannot accurately be accounted for at this stage.

Figure 6 shows the axial temperature profile of the CFD simulations compared to the measurement of the thermocouples in the external wall for case number 43. It is worth noting that thermocouples readings were shown for only 3 cases in [2]. As can be seen in the figure, the temperature distribution suggests the heat source is not uniform as it drops in the middle of the test section. However, it is not possible to recreate the experimental conditions given the lack of further details. Figure 6 shows the best combination of thermophysical properties along with the two combinations that represent the bounding cases. It can be seen that the uniform volumetric heat source approach is valid, especially with a lack of rationale for choosing a less standard approach to applying the thermal energy.

Due to the aforementioned spread in the wall temperature values, it was observed that the uncertainty in the material thermophysical properties could lead to a significant variance in the heat transfer properties of FLiNaK. This uncertainty coming from the material properties, added to the lack of reliable estimates of experimental uncertainties, impeded an appropriate validation of RANS models for molten salt flows. From this point onwards, it was assessed that a sensitivity analysis on how these properties’ individual variations affected heat transfer properties was needed.

Table 5

Important input parameters for the ten modelled cases.

thumbnail Fig. 5

Difference between the wall temperature measured in the experiment and the wall temperature resulting from the simulations for each of the ten cases. The various data points in each case correspond to different selections of material properties of FLiNaK.

thumbnail Fig. 6

Axial temperature on the outer wall of the test section measured with thermocuples in [2] and obtained in the simulations for case 43.

5 Sensitivity analysis

A widely used parameter to measure the sensitivity of a simulation result S(X) to changes in a parameter Xi is termed asthe sensitivity coefficient, whose definition is the following: (9)

where Xi is one element of the vector , which includes all dependencies of the variable S.

Using a Taylor series approximation, the uncertainty of a function S(X) with n uncorrelated parameters can be accounted for by using the following expression: (10)

For simple cases of an algebraic model, these sensitivity coefficients may be calculated analytically. However, the most common case is that the model is a complex numerical simulation where no algebraic model can be used. Using this analysis as an example, a change in one of the material properties may produce a change in the velocity profile of the fluid and also in the heat transfer properties between the liquid and the solid. Because of the existence of these mentioned coupled effects, the sensitivity coefficient was calculated with data from computational fluid dynamics simulations.

The procedure used in this section was to run the simulation with nominal values of the parameter vector . For the next simulations, perturbed values for the input parameter Xi are used. Finally, using a forward finite difference approximation in the parameter space, the sensitivity coefficient is calculated as follows [24]: (11)

For this particular case, the simulation result of interest S(X) was the heat transfer coefficient between the molten salt and the solid pipe. The vector of parameters consisted of the four thermophysical properties of FLiNaK.

It is worth pointing out that the One Factor At a Time Approach (OAT) [25] was used in the following analysis. In this approach, only one parameter changes its value between consecutive simulations, and so, in a deterministic model, the analyst can determine exactly what effect is caused by changing the parameter. In order to obtain good results, the model must be linear in the sense that Gaussian distributions are assumed between input and output values. In addition, there should not be a significant covariance between the input parameters of the sensitivity analysis. Even though it is known that the thermophysical properties are related to similar underlying molecular mechanisms [26] and that the model may have non-linear effects, the current analysis presents a first approach to sensitivity analysis in molten salt flows. Therefore, for the sake of simplification, the model response is assumed linear and that the covariance is not significant. Additional studies should be conducted in the future to validate this claim.

Taking into account every correlation presented in Table 1, an average temperature dependent function was calculated for each property. In order to do this, the temperature dependent expression for each property was evaluated at each temperature in the range of 800–1100K. Next, these values were averaged, yielding one average value of the property for each temperature. Finally, a fitting analysis was made to obtain the average function which better described the behaviour of the material property in question (see Fig. 7). From this point onwards, in each one of the simulations, only one of the properties was varied while keeping this same functional shape obtained in the fitting. The variations for the analysis were constant steps and the range of these variations was in agreement with respect to maximum and minimum values of the property of interest (as shown in Fig. 7). The remaining average temperature dependent properties were kept invariant for each analysis. As an example, the input viscosities for the sensitivity analysis are shown in Figure 8.

In Figure 7 the four models for viscosity in addition to the average viscosity in the temperature range of interest are plotted. From Figure 8, the average viscosity function and the simulation inputs for the sensitivity analysis canbe observed. An analogous process of obtaining an average function and adding constant steps to build the simulation inputs was done for the density, the thermal conductivity and the specific heat. Particularly for the thermal conductivity case Smirnov’s correlation [7] (shown in Tab. 1) was used, and the variations were done according to the studies done by Williams et al. [12], which states that the values should range between 0.6 W/m K and 1.0 W/m K at 973 K.

As it was mentioned before, the sensitivity analysis in this section was targeted at obtaining sensitivity coefficients for the heat transfer coefficient variable. Thus, each one of the eight black curves labeled as Simulation Inputs in Figure 8 were used inthe simulations in order to obtain the corresponding heat transfer coefficients as a function of the property value evaluated atthe bulk temperature. The simulations were done in three different temperature ranges, determined by different temperatureinlet boundary conditions. Hence, some of the typical temperatures of operation present in a MSR were covered.

The resulting heat transfer coefficients obtained in the simulations as a function of viscosity, density, specific heat and thermal conductivities are plotted in Figure 9. Taking as an example the viscosity sensitivty analysis, from Figure 9d it can be observed that there are eight different values of heat transfer coefficients for each temperature range. This matches the number of simulation inputs described in Figure 8. From the data shown in Figure 9, the results of the calculation of the sensitivity coefficients for different temperature ranges are shown in Table 6.

Note that the sensitivity coefficient is positive for three of the four properties considered with viscosity alone having a negative value. This is not surprising since an increase in viscosity would decrease the Reynolds number, thereby diminishing the amount of turbulence as well as the transfer of thermal energy.

Nusselt number calculations using Dittus Boelter and Sieder Tate correlations were performed in order to compare and assess the validity of using both of these correlations. The comparison was done using the results for the viscosity sensitivity analysis, which was arepresentative case (see Fig. 10).

As it can be seen from Figure 10, Sieder Tate predictions agree much better with CFD simulation results compared to Dittus Boelter. This is not surprising in that the Sieder Tate formula accounts for the variability in fluid properties near the wall, which can influence the heat transfer performance. The slope of the curves in Figure 10 also reveal a good agreement between Sieder Tate and our CFD simulations suggesting similar sensitivity from the viscosity on the heat transfer coefficient could be predicted based on the Sieder Tate correlation. It is important to note that the CFD simulations allow a direct calculation of the surface temperature, and therefore provide the evidence of the validity of the Sieder Tate formula in this scenario. Additional work should be performed with other turbulence models and higher fidelity simulations to more completely establish this validity.

Apart from the sensitivity coefficients presented above in equation (11), a dimensionless sensitivity parameter is an alternative way to provide information about how a certain parameter affects the variable of interest. Consequently, the dimensionless sensitivity coefficient is defined as: (12)

Equation (12) suggests a rather simple interpretation for S*, namely a fractional variation in the calculated observable quantity (in this case the heat transfer coefficient) that is produced by a specific fractional change in a selected input parameter (in this case material properties). In equation (12), the observable quantity h is defined as the heat transfer coefficient values obtained by using the average values of the material properties x. For the viscosity sensitivity analysis, the heat transfer coefficient h appearing in equation (12) would be calculated with the average viscosity red curve in Figure 8 and the remaining property average values. The amount of variation of x is δx while the resulting change in h is δh. Using these definitions, the dimensionless sensitivity coefficients were calculated and are presented in Table 7. Also, the maximum variations in percentage of a material property according to the different correlations can be observed.

According to the results shown in Table 7, viscosity, thermal conductivity and density have the highest dimensionless sensitivity coefficient. Nevertheless, the percentage change in density is very small compared to the thermal conductivityand viscosity, making it less relevant in terms of how it might affect the prediction of the heat transfer coefficient. This is due to the fact that the determination of the value of density is a less complex task than the measurement of the other two properties. In contrast, the thermal conductivity values reported by [6] and [12], may vary between 0.6 and 1.0 W/m K approximately, which accounts for a 26% spread.

Therefore, to be able to producing high quality experimental data for future validation purposes, a more confident understanding of material properties is fundamental. More specifically to future molten salt efforts, thermal conductivity and viscosity need to be accurately described since they are the most sensitive parameter.

thumbnail Fig. 7

Temperature dependent correlations for viscosity and the average viscosity function obtained by exponential fitting.

thumbnail Fig. 8

Viscosity inputs for the sensitivity analysis CFD simulations.

thumbnail Fig. 9

Results of the sensitivity analysis 34.

Table 6

Sensitivity coefficients for each property in different temperature ranges.

thumbnail Fig. 10

Heattransfer coefficient versus thermal conductivity values obtained from the sensitivity analysis. Comparisons with Dittus Boelter and Sieder Tate weremade.

Table 7

Dimensionless sensitivity coefficients for each property, and their respective maximum variation determined by the different correlations.

6 Conclusions

Reynolds averaged Navier Stokes using the kω SST model was used to perform computational fluid dynamics simulations attempting to reproduce the experiment results obtained by Grele and Gedeon [2] in 1954. Lack of detail and uncertainty estimates of the experimental data prevented a full validation effort. Despite the inability to conducting full validation studies, the influence of the variation of the thermophysical properties was gauged as they pertained to external wall temperatures and heat transfer coefficients. A large spread of values was observed for the external wall temperature of the simulations.

In light of the uncertainty of FLiNaK’s thermophyisical properties, a sensitivity analysis using a one-factor at a time approach (OAT) was done to account for the influence of these variations on the heat transfer coefficient for each one of the properties. The analysis proved that viscosity and thermal conductivity are the most crucial properties for conjugate heat transfer simulations using FLiNaK. The simulations also demonstrated that Sieder Tate correlation can be used with a reasonable level of trust since the results agreed well with that expression. It is not surprising that the Dittus Boelter expression did not agree as well, especially considering it cannot account for different viscosity at the wall compared to the bulk. This considered analysis can be extended to other molten salts due to the similarity in the dimensionless numbers (Prandtl and Reynolds) used in the study.

This work provides evidence on the importance of conducting further research in the material properties for FLiNaK, as they have a significant effect in the prediction of the heat transfer coefficient. Additional turbulence models should also be run to further validate the use of the Sieder Tate correlation and gain further insight into the flow physics.

For future studies, a sensitivity analysis using Large Eddy Simulations (LES) and DNS is suggested with the aim of reducing the impact arising from the error in the turbulence model used and the assumptions built into those models. In addition, it is important to note that this entire analysis was based on main effects of the material properties. Further studies to quantify any interactions between the variations of more than one property at the same time are suggested. This potential analysis would serve to quantify the covariance between the uncertainty of two properties. In this way a probability density functions of inputs and outputs in the current sensitivity study could be characterized in order to obtain the non-linearity effects of the thermophysical properties’ effect on the heat transfer coefficient.

Author contribution statement

Ramiro Freile has written the article. Mark Kimber has contributed to this work by providing financial support through a federal grant and advising on the different topics discussed.

References

  1. H. Bussier, S. Delpech, V. Ghetta, D. Heuer, D.E. Holcomb, V. Ignat’ev, E. Merle-Lucotte, J. Serp, The molten salt reactor (MSR) in generation IV: overview and perspectives, Prog. Nucl. Energy 77, 308 (2014) [CrossRef] [Google Scholar]
  2. M.D. Grele, L. Gedeon, Forced-convection heat-transfer characteristics of molten FLiNaK flowing in an Inconel X system (National Advisory Committee for Aeronautics, 1954) [Google Scholar]
  3. H.W. Hoffman, J. Lones, Fused Salt Heat Transfer Part II: Forced Convection Heat Transfer in Circular Tubes Containing NaF-KF-LiF Eutectic, ORNL-1977, Oak Ridge National Laboratory, 1955 [Google Scholar]
  4. I.B. Vriesema, Aspects of Molten Fluorides as Heat Transfer Agents for Power Generation, WTHD No. 112, Delft University of Technology, 1979 [Google Scholar]
  5. V. Ignat’ev, et al., Heat Exchange During the Flow of a Melt of LiF-NaF-KF Fluoride Salts in a Circular Tube, Sov. At. Energy 57, 123 (1984) [Google Scholar]
  6. C.T. Ewing, et al., Radiant transfer of heat in molten inorganic compounds at high temperatures, J. Chem. Eng. Data 7,246 (1962) [CrossRef] [Google Scholar]
  7. M.V. Smirnov, V.A. Khoklov, E.S. Filatov, Thermal conductivity of molten alkali halides and their mixtures, Electrochim. Acta 32, 1019 (1987) [CrossRef] [Google Scholar]
  8. J. Ambrosek, M. Anderson, K. Sridharan, T. Allen, Current status of knowledge of the fluoride salt (FLiNaK) heat transfer, Nucl. Technol. 165, 166 (2009) [CrossRef] [Google Scholar]
  9. M.S. Sohal, M.A. Ebner, P. Sabharwall, P. Sharpe, Engineering database of liquid salt thermophysical and thermochemical properties (No. INL/EXT-10-18297), Idaho National Laboratory (INL), 2010 [CrossRef] [Google Scholar]
  10. C. Xu, Z. Wang, Y. He, X. Li, F. Bai, Sensitivity analysis of the numerical study on the thermal performance of a packed-bed molten salt thermocline thermal storage system, Appl. Energy 92, 65 (2012) [CrossRef] [Google Scholar]
  11. P. Sabharwall, E.S. Kim, M. McKellar, N. Anderson, Process heat exchanger options for the advanced high temperature reactor (No. INL/EXT-11-21584), Idaho National Laboratory (INL), 2011 [Google Scholar]
  12. D.F. Williams, L.M. Toth, K.T. Clarno, Assessment of Candidate Molten Salt Coolants for the Advanced High-Temperature Reactor (AHTR), ORNL/TM-2006/12, Oak Ridge National Laboratory, Oak Ridge, TN, 2006 [CrossRef] [Google Scholar]
  13. G.J. Janz, R.P.T. Tomkins, Physical Properties Data Compilations Relevant to Energy Storage: IV Molten Salts: Data on Additional Single and Multi-Component Salt Systems, National Standard Reference Data System, National Bureau of Standards Report NSRDS-NBS 61 Part IV, 1981 [Google Scholar]
  14. T. Allen, Molten Salt Database, Nuclear Engineering and Engineering Physics Department, University of Wisconsin (2010), http://allen.neep.wisc.edu/shell/index.php/salts [Google Scholar]
  15. F. Menter, Zonal two equation kw turbulence models for aerodynamic flows, in 23rd fluid dynamics, plasmadynamics, and lasers conference (1993), p. 2906 [Google Scholar]
  16. W.S. Kim, et al., Performance of a variety of low Reynolds number turbulence models applied to mixed convection heat transfer to air flowing upwards in a vertical tube, Proc. Inst. Mech. Eng. C 218, 1361 (2004) [CrossRef] [Google Scholar]
  17. F. Menter, J. Carregal Ferreira, T. Esch, B. Konno, The SST Turbulence Model with Improved Wall Treatment for Heat Transfer Predictions in Gas Turbines, in Proceedings of the International Gas Turbine Congress, Tokyo, 2003 [Google Scholar]
  18. D.C. Wilcox, in Turbulence modeling for CFD (DCW industries, La Canada, CA, 1998), Vol. 2, pp. 172–180 [Google Scholar]
  19. Y. Chen, Z. Tang, N. Wang, Numerical prediction of turbulent convective heat transfer with molten salt in circular pipe, in NURETH-16, Chicago, IL, August 30-September 4, 2015 [Google Scholar]
  20. Y.M. Ferng, K.-Y. Lin, C.-W. Chi, CFD investigating thermal-hydraulic characteristics of FLiNaK salt as a heat exchange fluid, Appl. Thermal Eng. 37, 235 (2012) [CrossRef] [Google Scholar]
  21. F. Moukalled, L. Mangani, M. Darwish, An Advanced Introduction with OpenFOAM and Matlab, in The finite volume method in computational fluid dynamics (2016), pp. 3–8 [CrossRef] [Google Scholar]
  22. E.S. Chaleff, T. Blue, P. Sabharwall, Radiation heat transfer in the Molten Salt FLiNaK, Nucl. Technol. 196, 53 (2016) [CrossRef] [Google Scholar]
  23. Inconel X-750 Technical Data. High Temp Metals, 2015, www.hightempmetals.com/techdata/hitempInconelX 750data.php. Accessed November 2018 [Google Scholar]
  24. American Society of Mechanical Engineers, Standard for Verification and Validation in Computational Fluid Dynamics and Heat Transfer: An American National Standard, American Society of Mechanical Engineers, 2009 [Google Scholar]
  25. A. Saltelli et al., Global Sensitivity Analysis (John Wiley Sons, 2008) [Google Scholar]
  26. A.E. Gheribi, D. Corradini, L. Dewan, P. Chartrand, C. Simon, P.A. Madden, M. Salanne, Prediction of the thermophysical properties of molten salt fast reactor fuel from first-principles, Mol. Phys. 112, 1305 (2014) [CrossRef] [Google Scholar]

Cite this article as: Ramiro Freile and Mark Kimber. Influence of molten salt-(FLiNaK) thermophysical properties on a heated tube using CFD RANS turbulence modeling of an experimental testbed, EPJ Nuclear Sci. Technol. 5, 16 (2019)

All Tables

Table 1

Correlations for thermophysical properties of FLiNaK taken from [9].

Table 2

Hydrodynamic boundary conditions used.

Table 3

Energy boundary conditions used for the fluid region.

Table 4

Energy boundary conditions used for the solid region.

Table 5

Important input parameters for the ten modelled cases.

Table 6

Sensitivity coefficients for each property in different temperature ranges.

Table 7

Dimensionless sensitivity coefficients for each property, and their respective maximum variation determined by the different correlations.

All Figures

thumbnail Fig. 1

Experimental results obtained by Grele and Gedeon [2] and Hoffman and Loans [3] compared with Dittus Boelter correlation.

In the text
thumbnail Fig. 2

Experimental results obtained by Vriesema [4] compared with Dittus Boelter.

In the text
thumbnail Fig. 3

Grele and Gedeon, Hoffman and Loans and Vriesema results for Inconel piping reanalyzed by Ambrosek et al. [8].

In the text
thumbnail Fig. 4

Perpendicular section of the computational domain used in the simulations.

In the text
thumbnail Fig. 5

Difference between the wall temperature measured in the experiment and the wall temperature resulting from the simulations for each of the ten cases. The various data points in each case correspond to different selections of material properties of FLiNaK.

In the text
thumbnail Fig. 6

Axial temperature on the outer wall of the test section measured with thermocuples in [2] and obtained in the simulations for case 43.

In the text
thumbnail Fig. 7

Temperature dependent correlations for viscosity and the average viscosity function obtained by exponential fitting.

In the text
thumbnail Fig. 8

Viscosity inputs for the sensitivity analysis CFD simulations.

In the text
thumbnail Fig. 9

Results of the sensitivity analysis 34.

In the text
thumbnail Fig. 10

Heattransfer coefficient versus thermal conductivity values obtained from the sensitivity analysis. Comparisons with Dittus Boelter and Sieder Tate weremade.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.