Numerical investigation on heat transfer characterization of liquid lithium metal in pipe

Nomenclature
Cp	specific heat at constant pressure (J kg−1·k−1)
h	convective heat transfer coefficient (W m−2·K−1)
k	turbulent kinetic energy
L	pipe length (m)
Nu	Nusselt number
Pr	Prandtl number
P	pressure (Pa)
Pe	Peclet number
q	heat flux (W m−2)
Re	Reynolds number
r	pipe radius (m)
T	temperature (K)
$\vec v$	velocity (m s−1)
Greek
υ	momentum diffusion coefficient
α	diffusion coefficient
ω	dissipation rate
μ	dynamic viscosity (Pa s)
λ	thermal conductivity (W m−1·K−1)
ρ	density (kg m−3)
Subscripts
M	molecule
In	inlet
t	turbulence
W	wall
x	local position of pipeline

 | Show Table

DownLoad: CSV

1.   Introduction

The power supply of space nuclear reactors is expected to satisfy the requirements of manned space flight and deep space exploration, and liquid metals are primarily used in blanket nuclear fusion reactions as it can accommodate the extreme thermal loads (MW m⁻²) involved in the operation of space nuclear reactors^[1-3]. For example, lead (Pb) or bismuth–lead eutectic (LBE) are used in accelerator-driven subcritical systems and lead-alloy cooled fast reactors because of their high atomic number, favorable neutron properties, low melting point, and chemical stability^{[4, 5]}. Liquid NaK has a wide liquid temperature range (262~1073 K) and low vapor pressure in a vacuum environment^[6]. Metallic sodium (Na) with an extremely high thermal conductivity has been used in the design of Phénix reactors^[7]. Liquid Li is planned for use in sp-100 space reactor power supplies because of its high boiling temperature^[8].

In recent decades, the adaptation of liquid Na and LBE to nuclear reactors has been investigated extensively^[9-11]. In terms of Pb or LBE, Ma et al.^[12] investigated the hydrodynamic performance of liquid LBE flowing through a straight tube and U-type heat exchangers under a low Reynolds number, from which the pressure drop, friction coefficient, heat transfer coefficient, and other experimental data were obtained. In addition, Wang et al.^[13] numerically investigated the effect of buoyancy force on the convective heat transfer of LBE in a three-dimensional circular tube at different tilt angles and aspect ratios, and the results showed that a smaller angle of inclination and larger diameter resulted in better convective heat transfer. For liquid Na or NaK, gas entrainment in the circuit is vital to the safety of nuclear reactors^[14]. Xiao et al.^[15] experimentally investigated the heat transfer of liquid Na at a high temperature (300～700 ℃) and a low Peclet(Pe) number (20～70), as well as at a low temperature (250～270 ℃) and a high Pe (125～860). Additionally, theoretical and experimental analyses pertaining to the superheating of liquid Na in annulus boiling were performed, and the obtained semi-empirical formulas agreed well with the experimental data. Wang et al.^[16] simulated the local turbulent heat transfer of liquid Na in a circular tube, including the effects of turbulence, heating, and geometry. They discovered that the inlet velocity and the gap size of the annulus significantly affected the turbulent heat transfer of liquid Na, whereas the inlet temperature and the density of the heated flow imposed limited effects.

Compared with Pb, LBE, and Na, liquid Li presents higher adaptability to temperatures up to 1600 K, lower density, and more favorable neutron characteristics; therefore, it can satisfy the cooling requirements of megawatt-class space nuclear reactor power supplies. Xu et al.^[17] presented a thermoelectric magnetohydrodynamic driven flowing liquid Li in an open surface stainless steel trench structure,which can withstand heat fluxes of up to 10 MW·m⁻² for 10 s without significant evaporation. However, because high-temperature liquid metals present a high risk during experiments, the measurement of their physical parameters such as speed, temperature, and pressure is extremely difficult^[18]. Hence, numerical simulation is an effective option for investigating the heat transfer characteristics of flowing liquid metals. Ziyaddin et al.^[19] performed a numerical analysis to investigate the effect of magnetic field on the laminar convective heat transfer of liquid Li. It was discovered that as the magnetic field increased, the heat transfer performance enhanced. However, as the physical properties of liquid Li change significantly with temperature, the thermal properties of liquid Li are regarded as constant physical properties in most studies, and the flow of liquid Li in practical engineering is typically turbulent. When considering the characteristics of a space environment, whether knowledge regarding relevant liquid metals can be used for Li cooling requires further verification. It is reasonable to expect that the numerical analysis of the turbulent heat transfer of liquid Li can facilitate its application in space nuclear reactors. However, to the best of our knowledge, the turbulent flow of liquid Li in a space environment has not been investigated experimentally or numerically.

Herein, a steady-state, two-dimensional axisymmetric model is used to describe the flow and heat transfer processes of liquid Li in a straight tube used in the primary circuit of space nuclear reactors. First, the heat transfer properties of liquid Li, Na, and LBE were compared at different inlet temperatures. Subsequently, the effects of different parameters, e.g., the flow rate, inlet temperature, and heat flux applied to the pipe wall, on the heat transfer of liquid Li are presented in detail. This study provides a basis for the flow and heat transfer characterization of high-temperature liquid Li, reveals the mechanism of internal heat flux, and facilitates the advancement of Li-cooled nuclear reactors.

2.   Model development

Fig. 1 illustrates a liquid metal flowing in a pipe with a radius of 0.025 m and length of 3 m. A constant heat flux was imposed around the pipe to heat the liquid metal. In the following calculations, one-half of the pipe was selected as the flow calculation area, and the intermediate symmetry axis was used as the 0-axis to analyze the results.

Figure  1.  Schematic illustration of liquid metal flowing in a pipe.

DownLoad: Full-Size Img PowerPoint

To simplify the calculations, the following were assumed:

(Ⅰ) The thickness of the pipe wall is negligible.

(Ⅱ) The boundary is non-slip.

(Ⅲ) The liquid metals are incompressible fluids.

(Ⅳ) The flow of the liquid metal is a single-phase flow without phase transition.

(Ⅴ) The effects of gravity and buoyance are negligible (since a space environment is assumed).

2.1   Governing equations

The standard k–ω model was used to simulate the turbulent flow of the liquid metal in a pipe^[20].The relevant control equations for the flowing liquid metal are as follows:

The continuous equation is expressed as

where ρ and $\vec v$ denote the density and flow velocity of the liquid metal, respectively.

The momentum equation is expressed as

where P denotes pressure, μ viscosity, and $\vec F$ volume force.

The turbulent kinetic energy k and kinetic energy dissipation rate ω are expressed as follows^[21]:

where μ_T is the turbulent viscosity, expressed as ${\mu _T} = \rho \dfrac{k}{\omega }$; σ_k^*, σ_ω, ξ, β₀^*, and β_o are the system turbulence model constants.

The energy equation is expressed as

where C_p is the specific heat capacity at a constant pressure, λ is the thermal conductivity, and Q is the heat source.

2.2   Relevant parameters

The initial turbulent kinetic energy k₀ and initial dissipation rate ω₀ can be expressed as follows:

where I_T and L_T are the initial turbulence intensity and initial turbulence length, respectively.

The convective heat transfer coefficient h_x can be expressed as

where T_w,x and T_b,x indicate the wall temperature and liquid metal bulk temperature at position x along the heat pipe, respectively.

The dimensionless Nusselt (Nu) number, Prandtl (Pr) number,and Pe can be expressed as follows:

where υ is the momentum diffusion coefficient and α is the thermal diffusion coefficient.

The radial turbulent heat flux q_t is expressed as

where v_r′ is the fluctuating parts of the radial velocity field,and T′ is the fluctuating region of the temperature field; α_t is the thermal diffusion coefficient of turbulence.

The molecular heat flux q_M is expressed as

where α_M is the thermal diffusion coefficient of the molecule. Hence, the coefficient of turbulence thermal diffusion to molecular thermal diffusion can be expressed as

where q_w=q_t+q_M, and q_w is the constant wall heat flux.

2.3   Boundary conditions and turbulent Pr model

The boundary conditions of the velocity inlet and pressure outlet were adopted. The symmetry condition was applied on two symmetry planes. Because the k–ω model is a low Reynolds number model, the dimensionless distance Y⁺ of the first grid must be verified. The Y⁺ of the first layer grid was between 0.92 and 0.98, which satisfies Y⁺≤1. Hence, the default values were assigned to the inlet turbulent kinetic energy and specific dissipation rate, and a constant heat flux was applied to the wall.

Liquid metals typically have lower Pr values than general fluids. To improve the accuracy of the simulation, the default turbulent planetary Prandtl number Pr_t of the model must be corrected. Cheng and Tak^[22] provided a modified formula to obtain the Pr_t of LBE, as follows:

In the numerical simulation of liquid Li and Na, the equations above were applied to obtain the corrected Pr_t (i.e., instead of using the default Pr_t). To simplify the calculation, Pe was determined from the inlet temperature to obtain the modified Pr_t.

3.   Numerical results

3.1   Calculation details

All physical parameters in the calculations were provided as a function of temperature. The main physical parameters, including density, constant pressure specific heat, dynamic viscosity, and thermal conductivity, are listed in Table 1.

Table  1.  Physical parameters of Li, Na, and LBE^{[23, 24]}.

Metal	Properties	Expression
Li(454～1608 K)	ρ /(kg m−3)	454～900 K: ρ=562−0.1·T
		900～1608 K: lnρ=6.0669+8.3105×10−4·T−1.5042×10−6·${T^2}$+ 1.0686×10−9·${T^3}$−3.7709×10−13·${T^4}$+5.1832×10−17·${T^5}$
	Cp / (J kg−1·k−1)	454～900 K: Cp=0.0014·${T^2}$−2.4066·T+5183.4
		900～1608 K: lnCp=1000×(1.4122×10−1+5.1039×10−3·T−7.8268×10−6·${T^2}$+5.7796×10−9·${T^3}$− 2.0587×10−12·${T^4}$+2.8515×10−16·${T^5}$)
	μ / (Pa s)	454～900 K: $\log \mu = - 3.080 + \dfrac{ {57.63} }{T} - 5.172 \times {10^{ - 4} } \cdot T$
		900～1608 K:ln(104·μ)=9.4544−2.7505×10−2·T+ 3.6732×10−5·${T^2}$−2.5268×10−8·${T^3}$+ 8.6852×10−12·${T^4}$−1.1841×10−15·${T^5}$
	λ /(W m−1·K−1)	454～900 K: λ=21.874+0.056225·T−1.8325×10−5·${T^2}$
		900～1608 K: lnλ=4.2556−2.5537×10−3·T+5.0480×10−6·${T^2}$−4.0821×10−9·${T^3}$+1.5784×10−12·${T^4}$−2.3730×10−16·${T^5}$
Na(375～1255 K)	ρ /(kg m−3)	$\rho = 219 + 275.32 \times \left( {1 - \dfrac{T}{ {2503.7} } } \right) + 511.58 \times {\left( {1 - \dfrac{T}{ {2503.7} }} \right)^{0.5} }$
	Cp /(J kg−1·k−1)	Cp=1658.2−0.84790·T+4.4541×10−4·${T^2}$−2.9926×10−6·${T^{ - 2}}$
	μ /(Pa s)	$\ln \mu = - 6.4406 - 0.3958 \times \ln T + \dfrac{ {556.835} }{T}$
	λ /(W m−1·K−1)	λ=124.67−0.11381·T+5.5226×10−5·${T^2}$−1.1842×10−8·${T^3}$
LBE(400～1100 K)	ρ /(kg m−3)	ρ=11096−1.3236·T
	Cp /(J kg−1·K−1)	Cp=159−2.72×10−2·T+7.12×10−6·${T^2}$
	μ /(Pa s)	λ=3.61+1.517×10−2·T+7.12×10−6·${T^2}$
	λ /(W m−1·K−1)	$\mu = 4.94 \times {10^{ - 2} }\cdot\exp \left( {\dfrac{ {754.1} }{T} } \right)$

 | Show Table

DownLoad: CSV

The calculations were performed based on the finite element method using COMSOL multiphysics. A non-isothermal heat transfer multiphysics field was constructed using a turbulence module and a fluid heat transfer module. The grid contained 128000 elements, and the calculated residuals were set to 10⁻³.

3.2   Model validation

The accuracy of the model was verified by selecting LBE as the working fluid and comparing the simulated Nu with the Lyon empirical correlation as well as the experimental data by Johnson et al.^[25], Lyon correlated Nu as a function of Pe as

The numerical simulation results of Mochizuki were similar to those obtained using the Lyon empirical correlation and those reported in Ref. [26].

The variation of Nu with Pe is illustrated in Fig.2. As shown, for all the Pe numbers investigated, the numerically simulated Nu values were similar to those from the Lyon empirical correlation and the experimental data of Johnson et al.. This indicates that the model can be employed in the present study to simulate liquid metals.

Figure  2.  Comparison of simulation Nu, obtained from simulation, Lyon empirical correlation and experimental data of Johnson et al.

DownLoad: Full-Size Img PowerPoint

4.   Results and discussion

To compare the heat transfer characteristics of liquid Na, LBE, and Li, the coefficients of local convective heat transfer of these three liquid metals on the wall surface were investigated using the same inlet temperatures. Fig. 3 a shows that at 600 K, the convective heat transfer coefficients of liquid Na and Li, did not differ significantly, whereas that of the LBE was the smallest value. Therefore, compared with liquid Na and LBE, liquid Li did not offer significant advantages in terms of heat transfer at lower inlet temperature ranges. As shown in Fig. 3 b, when the inlet temperature increased to 1000 K, the local convective heat transfer coefficients of the liquid LBE and Li increased. For Li, the value increased from 9702 to 10200 W m⁻²·K⁻¹. However, the convective heat transfer coefficient of liquid Na decreased from 10000 to 7200 W m⁻²·K⁻¹. Therefore, at T_in = 1000 K, liquid Li exhibited the best heat transfer performance among the three liquid metals, demonstrating its advantages in high-temperature applications. In the following sections, we focus on the heat transfer characteristics of liquid Li.

Figure  3.  Heat transfer coefficients of three types of liquid metals at inlet temperatures of (a) 600 K and (b) 1000 K.

DownLoad: Full-Size Img PowerPoint

4.1   Effects of inlet velocity

The effect of the inlet velocity on the heat transfer of liquid Li, was investigated based on inlet velocities of 0.2, 0.4, and 0.6 m/s. The inlet temperature was set to 600 K, and the wall heat flux was set to 400000 W m⁻²·K⁻¹. As shown in Fig. 4 a, the liquid Li reached the full development stage after passing through the inlet of the pipeline for a certain distance, and its convective heat transfer coefficient stabilized gradually after passing through the inlet area. As the inlet velocity increased, the convective heat transfer coefficient at the outlet increased correspondingly, reaching 9704, 12099, and 14452 W m⁻² k⁻¹ at velocities of 0.2, 0.4, and 0.6 m/s, respectively.

Figure  4.  Effects of inlet velocity on heat transfer: (a) Local convective heat transfer coefficients along pipeline. (b) Turbulent heat flux. (c) Ratio of turbulent diffusion coefficient to molecular diffusion coefficient along the radius of the pipe.

DownLoad: Full-Size Img PowerPoint

Heat flux affects the numerical simulations of liquid metal flow and heat transfer in pipes. As shown in Fig. 4 b, in the radial direction of the outlet pipe, the turbulent heat flux increased from the center to the wall of the pipe first, and then decreased gradually after reaching a peak value, reflecting a parabolic distribution. As the inlet velocity increased, the peak value of the radial turbulent heat flux increased gradually, and the position gradually approached the pipe wall. This is because the increase in velocity resulted in increased turbulent disturbance near the wall, which consequently intensified turbulent heat transfer.

As shown in Fig. 4 c, the distribution of α_t/α_M from the center of the pipe to the wall was similar to the turbulent heat flux distribution, which presents a parabolic distribution along the radius of the outlet, and the peak value increased with the velocity. Owing to the low Pr of the liquid Li, at the inlet velocities of 0.2 and 0.4 m s⁻¹, the value of α_t/α_M was less than 1.0, which implies that turbulent heat transfer was less significant than molecular heat transfer. Hence, molecular heat transfer dominated the radial heat transfer of liquid Li. When the inlet velocity increased to 0.6 m/s, the peak value of α_t/α_M exceeded 1.0, which implies that turbulent heat transfer dominated in the radial direction of the pipe, whereas molecular heat transfer dominated in the other directions.

4.2   Effects of inlet temperature

The thermal properties of liquid Li changed with temperature, as indicated in Table 1. Hence, temperature is vital to the heat transfer of liquid Li. Fig. 5 a shows that the convective heat transfer coefficients were approximately 14000 W m⁻²·K⁻¹ when the inlet temperature increased from 600 to 1000 K. At the inlet temperature of 1400 K, the convective heat transfer coefficient increased to 14960 W m⁻²·K⁻¹. Therefore, heat transfer was more significant at inlet temperatures exceeding 1000 K.

Figure  5.  Effects of inlet temperature on heat transfer performance: (a) Local convective heat transfer coefficients. (b) Ratio of turbulent thermal diffusion coefficient to molecular diffusion coefficient.

DownLoad: Full-Size Img PowerPoint

The effect of temperature on the ratio of turbulent to molecular diffusion coefficients is shown in Fig. 5 b. When the inlet temperatures were 1000 and 1400 K, the value of α_t/α_M was less than 1.0, indicating that molecular heat transfer dominated the radial heat transfer process. As the inlet temperature decreased to 600 K, the value of α_t/α_M exceeded those at 1000 and 1400 K, and the peak value exceeded 1.0, indicating that the radial molecular heat transfer of liquid Li in the straight tube weakened as the inlet temperature decreased.

4.3   Effects of inlet wall heat flux

Three different heat flux densities, i.e., 2×10⁵, 4×10⁵, and 6×10⁵ W m⁻² were applied to the wall. The inlet velocity of liquid Li and the inlet temperature were maintained at 0.6 m s⁻¹ and 600 K, respectively. As shown in Fig. 6 a, the convective heat transfer coefficients along the pipeline were 14306, 14452, and 14580 W m⁻²·K⁻¹. Therefore, the wall heat flux density affected the thermal physical properties of liquid Li, which consequently affected the convective heat transfer coefficient. The convective heat transfer coefficient increased with the wall heat flux density.

Figure  6.  Effect of heat flux on heat transfer performance: (a) Local convective heat transfer coefficient. (b) Ratio of turbulent thermal diffusion coefficient to molecular diffusion coefficient.

DownLoad: Full-Size Img PowerPoint

The changes in α_t/α_M at different heat fluxes are shown in Fig. 6 b. As the heat flux on the wall increased, the peak position did not change significantly, and the value near the peak position exceeded 1.0, indicating that turbulent heat transfer was dominant at this position, whereas molecular heat transfer was dominant at other positions. As the heat flux increased, the value of α_t/α_M decreased, which is the same effect as increasing the inlet temperature of the liquid metal; this occurred because the temperature of the liquid metal increased with the heat flux.

The analysis above indicates that the ratio of turbulent to molecular diffusion coefficients depends on the inlet velocity, inlet temperature, and wall heat flux, and it can be greater or less than 1.0. This conclusion contradicts that in Ref. [27]. Previous studies have suggested that turbulent heat transfer dominates in water. Conversely, in liquid metals, molecular heat conduction dominates instead of turbulent heat transfer. The conclusions above were not based on comparing water and liquid metals under the same conditions. Furthermore, the effects of the inlet temperature, inlet velocity, and wall heat flux were not considered. Based on the current numerical simulation results, increasing the inlet velocity, decreasing the inlet temperature, and decreasing the wall heat flux can effectively weaken the role of molecular heat transfer but enhance turbulent heat transfer owing to the low Pr of liquid Li.

5.   Conclusions

In this study, a two-dimensional, steady-state, and non-isothermal heat transfer model for liquid Li turbulent flow in a straight pipe was developed. After performing validation based on LBE, a comparison between the simulation results and results obtained using the empirical formula indicatedthat the model can accurately capture the heat transfer characteristics of liquid Li. Based on this model, a comprehensive analysis of the key parameters, i.e., the inlet velocity, inlet temperature, and wall heat flux, was performed. The results showed that the heat transfer of liquid Li was better at high temperatures. When the inlet temperature exceeded 1000 K, the local convective heat transfer coefficient of liquid Li was significantly higher those that of liquid Na and LBE by more than 38%. The heat transfer performance of liquid Li can be further improved by increasing the inlet velocity, inlet temperature, and wall heat flux. For the radial heat flux, the ratio of turbulent to molecular diffusion coefficients presented a parabolic-like distribution along the radius of the pipe. The peak value was governed by the inlet velocity, inlet temperature, and wall heat flux, which can be either greater or less than 1. Increasing the inlet velocity, decreasing the inlet temperature, and decreasing the wall heat flux can effectively weaken the role of molecular heat transfer but enhance turbulent heat transfer owing to the low Pr of liquid Li. In future studies, the numerical simulation of liquid Li using different calculation models should be conducted to improve the computational accuracy, flow, and heat transfer of liquid Li in pipes of different shapes, and bending angles should be additionally considered for engineering applications. This study provides a basis for the heat transfer characterization of high-temperature liquid Li, reveals the mechanism of internal heat flux, and facilitates the advancement of Li-cooled space nuclear reactors.sents a parabolic-like distribution along the radius of the pipe.The peak value depends on the inlet velocity, inlet temperature, and wall heat flux, which can be either greater than one or not. Increasing the inlet velocity, decreasing the inlet temperature, and decreasing the wall heat flux can effectively weaken the role of the molecular heat transfer but enhance the turbulent heat transfer due to the low Pr of liquid Li. Furthermore, future research work should focus on the numerical simulation of liquid Li with different calculation models to improve the computational accuracy, flow and heat transfer of liquid Li in pipes with different shapes, and bending angles should also be considered for engineering application. Thus, this work can offer the flow and heat transfer characterization of the high-temperature liquid Li, reveal the mechanism of internal heat flux, and facilitate the advancement of Li-cooled space nuclear reactors.

Acknowledgments

The work is supported by the National Key R&D Program of China(2018YFB1900602).