MODELING OF THE FLOW AND HEAT TRANSFER OF SUPERCRITICAL CO2 FLOWING IN SERPENTINE TUBES

As a non-flammable, non-toxic refrigerant, supercritical CO2 (ScCO2) has been increasingly used for heat transfer applications. In this study, the ScCO2 flow and heat transfer in a set of full-size three-dimensional serpentine tubes were modeled with different inner diameters and tube pitches. The standard k-epsilon model was used for the turbulence modeling. The results show the effect of the different tube inner diameters and tube pitches on the flow and heat transfer of ScCO2 for a given flow flux or inlet Reynolds number. The heat transfer coefficient decreases as both the tube pitch and the inner diameter increase for a given mass flow rate. However, for a given inlet Reynolds number, the heat transfer coefficient first increases but then decreases with increasing tube inner diameter. The effect of the flow direction on the heat transfer performance was also studied for various inlet conditions. Downward flow results in a higher heat transfer coefficient than upward flow for inner diameters larger than 0.5 mm and the buoyancy effect can be ignored for Bo1×10-7 for the conditions studied here. These results can be used to optimize the tube shape and size in heat exchanger designs.


INTRODUCTION
CO2 is a non-flammable, non-toxic refrigerant with many attractive properties such as its low global warming potential and zero ozone depletion potential. Hence, CO2 is a promising substitute for hydrofluorocarbons (HFCs) to address steadily increasing international environmental concerns (Lorentzen, 1994;Pearson, 2005).
Due to its low critical temperature (31.3°C), CO2 is normally used as in trans-critical cycles for heat pumps and automotive air-conditioning systems (Groll and Jun, 2007;Zhang et al., 2007) with the heat rejection taking place in a gas cooler in the pseudocritical region. The thermophysical properties of supercritical fluids show drastic variations near the pseudo-critical point resulting in significantly different flow and heat transfer characteristics from conventional single phase flow (Kim et al., 2004). Therefore, the heat transfer characteristics of supercritical CO2 (ScCO2) flows in tubes must be carefully studied to improve system efficiencies.
The main factors (physical properties, buoyancy, flow acceleration, inlet conditions, flow direction, and tube diameter) affecting the turbulent flow and heat transfer of ScCO2 in tubes have been extensively studied, especially for straight round tubes (Kuang et al., 2004;Liu et al., 2017;Shiralkar and Peter, 1970;Zhang et al., 2016). The results from these studies have shown that turbulent heat transfer to ScCO2 is very complex because of the extreme property variations. Many studies have indicated that serious heat transfer deterioration can occur at supercritical pressures which significantly increases the wall temperature. During heating of ScCO2 in tubes, the heat transfer coefficient can significantly decrease at high heat fluxes when the bulk temperature is below the critical temperature and the wall temperature is above the critical temperature (Shitsman, 1963). The deterioration is most significant when the ScCO2 * Institute of Thermal Science and Technology, Shandong University, Jinan, Shandong, 250061, China † Corresponding author. Email: lixf@sdu.edu.cn ‡ Corresponding author. Email: cheng@sdu.edu.cn pressure is near the critical pressure (Shiralkar and Peter, 1970). During heating, the temperature can vary greatly in the radial direction with the fluid near the wall changing from a liquid-like state to a gas-like state. These changes significantly reduce the thermal conductivity neat the wall which reduces the heat transfer. This phenomenon is similar to the heat transfer deterioration during the transition from nucleate boiling to film boiling (Shitsman, 1963). During cooling of ScCO2, the heat transfer is enhanced when the wall temperature is less than the critical temperature and the fluid bulk temperature is greater than the critical temperature. This phenomenon is attributed to the lower temperature, liquid-like layer near the tube wall which has a higher thermal conductivity (Hiroaki et al., 1971;Krasnoshchekov et al. 1970;Shitsman, 1963). The heat transfer deterioration is also related to the reduced turbulence intensity in the fluid near the wall. The low turbulence intensity near the wall is attributed to reduced shear stresses due to buoyancy forces arising from the density difference between the high-density mainstream and the lowdensity near-wall flow (Hall and Jackson, 1969).
The buoyancy effect and flow acceleration caused by the large radial and axial density gradients have different effects on the heat transfer coefficients for ScCO2 flowing in vertical tubes with different inner diameters (Fewster et al., 2004;Jackson and Hall, 1979;Li et al., 2010;Mceligot et al., 2004). The heat transfer characteristics of ScCO2 in vertical tubes with different inner diameters have been measured in various studies (Jiang et al., 2004;Jiang et al., 2013;Li et al., 2010). These experimental results have shown that buoyancy significantly affects the heat transfer for inner diameters greater than 1 mm. The flow acceleration also significantly affects the heat transfer coefficient for flows in small diameter vertical tubes. A dimensionless number, Bo, has been used to evaluate the buoyancy effect on the heat transfer for flows in the vertical direction (Jackson and Hall, 1979). The results show that the buoyancy reduces the heat transfer for upward flows in vertical

Frontiers in Heat and Mass Transfer
Available at www.ThermalFluidsCentral.org circular channels for 5.6×10 -7 ≤Bo≤1.2×10 -6 with the heat transfer deterioration gradually reduced as Bo increases for 1.2×10 -6 ＜Bo≤8× 10 -6 . The heat transfer is enhanced for Bo＞8×10 -6 . In addition, the Grashof number (Gr) increases sharply with increasing tube diameter for flows in horizontal tubes, so the Grashof number is used to characterize the flow and heat transfer in horizontal tubes . The Nusselt number (Nu) deceases with decreasing tube diameter at constant wall temperature for inner diameters ranging from 0.5 to 2.16 mm (Liao and Zhao, 2002). Furthermore, the location and the extent of the deterioration are sensitive to the tube geometry (Shiralkar and Peter, 1970). Therefore, the inner diameter is a key factor affecting the heat transfer of ScCO2 in tubes. However, there are few studies on the influence of the tube geometry on ScCO2 flow and heat transfer in curved tubes.
In addition to straight round tubes, other shapes of tubes are often used for heat exchange , such as elliptical tubes (Yang et al., 2018) and serpentine tubes. Serpentine tubes are widely used in heat transfer applications. Unlike in straight tubes, the heat transfer coefficients in serpentine tubes are significantly affected by secondary flows (Wen et al., 2007;Wu et al., 2007;Zhao et al., 2011). The flow of ScCO2 in a serpentine tube is more unstable than in a straight tube. Since the flow instabilities reduce the heat transfer deterioration (Shiralkar and Peter, 1970), serpentine tubes can improve the heat transfer coefficients for ScCO2 flows. The secondary flows are due to a pair of counter rotating vortices whose axes are parallel to the main flow. The primary mechanism responsible for the vortices is the interaction of the centrifugal force and the pressure gradient over the tube cross section. The centrifugal force within the tube is related to the tube shape and size in addition to the flow parameters and fluid properties. Previous results have indicated that the secondary flows improve the heat transfer. However, there are few studies on the effect of the tube geometry (inner diameter and tube pitch) on the flow and heat transfer for flows in serpentine tubes. Therefore, the ScCO2 flow and heat transfer need to be further studied in various serpentine tube sizes.
The flow and heat transfer in serpentine tubes were modeled in this study to study the effects of the inner diameters tube pitch and flow direction on the flow and heat transfer coefficient to further understand the mechanisms driving the heat transfer enhancement for ScCO2 flow in various vertical serpentine tubes.

Physical Models and Boundary Conditions
The ScCO2 flows in steel serpentine tubes were assumed to be heated at a constant heat generation rate. As shown in Fig. 1, all the serpentine tubes consisted of one 50 mm straight adiabatic inlet section, one 88 mm long serpentine section, and one 50 mm straight adiabatic outlet section. The straight adiabatic inlet section was used to ensure that the fluid was in fully developed region as it entered the serpentine section. A series of serpentine tubes were studied with 4, 8, and 16 mm tube pitches (D) and 0.5, 1, 1.5 and 2 mm inner diameters (d). The serpentine tube models are referred to as "d_ D_ serpentine tube" for the inner diameter and tube pitch. For example, the 0.5_4_serpentine tube refers to the serpentine tube having an inner diameter of 0.5 mm and a tube pitch of 4 mm. The gravitational acceleration (g=9.81 m/s 2 ) was set along the positive or negative Z axes to model or upward flows.

Mesh Independence
ICEM 14.5 was used to form the structured hexahedron meshes. The Otype grid generation method was used to refine the mesh near the wall as shown in Fig. 2 (a). Then, the whole mesh was created by extending the mesh in the initial cross section along the axis. (Fig. 2 (b)). In this work, 12 serpentine tubes with different inner diameters and tube pitches were studied. The mesh independence was verified for each model. For example, three different meshes were considered with 890,000 elements (Case 1), 1,160,000 elements (Case 2), and 1,580,000 elements (Case 3) for the mesh independence verification of 0.5_4_serpentine tube. The outer wall temperatures of the 0.5_4_serpentine tube predicted using the three meshes are shown in Fig.  3. The temperature profiles along the periphery of the serpentine section (illustrated by the yellow line on the physical model in Fig. 1) almost coincide. The average relative error of the outer wall temperatures between Cases 2 and 3 is 0.0065% and the relative error of the heat transfer coefficient between Cases 2 and 3 is within 0.1%. Thus, the mesh with 1,160,000 elements was used for the calculations.

Numerical Solution
Fluent 14.5 was used for the three-dimensional numerical simulations. The pressure-based Navier-Stokes solution algorithm was used with the SIMPLE scheme for the velocity and pressure coupling. The convection terms in the energy, momentum and turbulence equations were all discretized using the second-order upwind method. The solution was considered to be converged when the residuals were less than 1×10 -8 for the continuity, momentum and energy equations. Three temperature monitors were set to monitor the inlet section, serpentine section and outlet section temperatures.
The ScCO2 flowed into the serpentine tubes with a pressure of 7.65 MPa and a temperature of 297.5 K and was then heated by a constant heat generation rate of 97.6 MW/m 3 . The mass flow inlet boundary condition was applied at the tube inlet with the pressure outlet boundary condition at the tube outlet. The inlet conditions were divided into two types with a constant mass flow rate (0.65 kg/h) or a constant Reynolds number (3592). The no-slip boundary conditions were used for all the walls. The thermal conditions on inner walls of the straight pipe sections and all outer walls were all adiabatic, and the thermal condition on inner walls of the serpentine sections were coupled. Meanwhile, the serpentine sections have an internal heat source with a constant heat generation rate. The cases considered in the present study are described in Table 1. All the ScCO2 thermophysical properties near the critical point at 7.65 MPa were taken from the NIST REFPROP database with the variations as shown in Fig. 4. Large numbers of points were selected from the thermophysical property profiles in the temperature range of 280 K~ 390 K to create correlations that were then input into Fluent by user defined functions. The thermophysical properties of CO2 are influenced dramatically by both of temperature and pressure. But the largest pressure drop in present paper is small enough to ignore (0.009 MPa). Therefore, the influence of pressure change on thermophysical properties was not considered in the correlation. The thermophysical properties are calculated as 35

Turbulence Model Selection
The flow was considered as turbulent since the inlet Reynolds number is greater than 2200 under all conditions in this paper. The extremely large thermophysical property variations near the critical point significantly affect the flow and heat transfer characteristics of supercritical fluids which differ from those for constant property fluids. In addition, the buoyancy and centrifugal force effects also complicate the flow and heat transfer processes in serpentine tubes compared to those in straight tubes. Therefore, a suitable turbulence model must be selected to ensure the credibility of the numerical solutions. The predictions using various turbulence models are compared with previous experimental data (Xu et al., 2015) in Fig. 5 using a three-dimensional model with d=0.953 mm and D=8.01 mm. Six turbulence models were evaluated, the k-kl-ω model, SST k-ω model, standard k-ε model, low Re k-ε model, RNG k-ε model and realizable k-ε model. The results in Fig. 5 show that the standard k-ε model, the RNG k-ε model and the realizable k-ε model best predicted the experimental data for the inner wall temperature with maximum errors less than 1.5%. The standard k-ε model was selected for the calculations.

Data Processing Method
The average convective heat transfer coefficient, h, was calculated as: Where Ainner is the inner wall area of the serpentine section and ∆tm is the logarithmic mean temperature difference: Where Tin,f and Tout,f are the average fluid temperatures at the inlet and outlet cross sections and Tin,w and Tout,w are the average wall temperatures at the inlet and outlet cross sections.
Bo was used to evaluate the buoyancy effects. Where g is the gravitational acceleration, β is the volume expansion coefficient, λf is the thermal conductivity, and ν is the kinematic viscosity.

Temperature and Velocity Distributions along the Flow Direction
The temperature and velocity distributions along the center plane of the 1_4_serpentine tube for upward ScCO2 flow are shown in Fig. 6. The maximum velocities occur along the inside of the curve with the minimum velocities along the outside of the curves just after the wall curvature reverses directions. This contributes to the secondary flow impingement. The minimum and maximum outer wall temperatures are then due to the different turbulence intensities. In addition, the centrifugal force causes the lower temperature, denser working fluid to flow towards the outside of the bend.

Effect of the Inner Diameter and Tube Pitch on The Temperature Distribution
The effects of the inner diameter and the tube pitch on the flow and heat transfer of ScCO2 were investigated for 12 serpentine tubes with four inner diameters and three tube pitches for constant mass flow rate and constant Reynolds number inlet conditions.

Constant mass flow rate
The outer wall temperature distribution along the flow direction in the various serpentine tubes are shown in Fig.  7 for upward ScCO2 flow. The outer wall temperature of the serpentine tubes increases rapidly as the heat transfer coefficient decreases with larger tube pitches when the inner diameter was held constant because the larger tube pitch reduces the centrifugal force and the turbulence intensity. The influence of the tube pitch on the centrifugal force and the turbulence intensity can be compared by comparing the velocity distributions shown in Fig. 8 across a typical cross-section in the middle of serpentine tubes with different tube pitches. Two pairs of vortices develop in the serpentine tube due to the reversed bends as shown by the velocity distribution in Fig. 8(a). However, with the larger tube pitch, the frequency of the direction changes is reduced, so there is only one pair of vortices in the cross section when D=16 and the secondary flow and turbulence intensities are lower. The outer wall temperature curves have periodic peaks that increase with increasing inner diameter and tube pitch as the flow and heat transfer of the working fluid in the serpentine tubes is affected by the buoyancy and centrifugal force. The secondary flow caused by centrifugal force will enhance heat transfer. As the flow periodically reverses direction in the bends, the centrifugal force on the fluid flowing in the serpentine tube periodically reverses. The direction of the buoyancy force component perpendicular to the axis at the same or opposite direction with the centrifugal force which then periodically strengthens and weakens the effect of the centrifugal force on the heat transfer. The centrifugal force decreases with increasing tube pitch, meanwhile the buoyancy increases with increasing inner diameter. Therefore, the influence of these periodic changes in the flow direction on the heat transfer coefficient increase with increasing pitch and inner diameter. The smaller inner diameters lead to more uniform temperature distributions since the radial thermal resistance increases with increasing inner diameter which increases the temperature difference between the fluid temperature near the wall and the mainstream temperature. Therefore, the mainstream temperature near the entrance is lower with larger inner diameters and the fluid temperature near the wall more quickly increases to the pseudo-critical temperature. The outer wall temperature then stabilizes earlier and then rises again quickly as shown in Fig. 7. Thus, a smaller inner diameter delays the rapid temperature rise and increases the temperature uniformity. With larger inner diameters, the rapid temperature rise occurs earlier and the temperatures are less uniform. Thus, the inner diameter, pitch and length of serpentine tubes should be selected to maximize the heat transfer and reduce the heat transfer deterioration.

Constant inlet Reynolds number The effects of the inner diameter
and tube pitch on the flow and heat transfer were also investigated with the inlet Reynolds number held constant. The outer wall temperature distribution and the centerline CO2 temperature distribution are shown in Fig. 9. The outer wall and centerline fluid temperatures both increase as the inner diameter decrease because the smaller inner diameter gives a smaller heat transfer area and a large heat flux since the total heat load is held constant. In addition, the mass flux decreases as the inner diameter is reduced since the Reynolds number is constant. Therefore, the fluid is heated more which leads to higher outer wall surface temperatures. In addition, the effect of the increasing inner diameter on the temperature distribution decreases with increasing inner diameter.  The temperature distribution for the 0.5 mm inner diameter serpentine can be divided into three stages while the temperature distributions in the other serpentine tubes only have two stages. The three stages are due to the dramatic changes in the specific heat (cp) of CO2 near the pseudo-critical point. Near the inlet, the fluid temperature near the wall is lower than the quasi-critical temperature and cp is small near the wall. The fluid then only absorbs a small amount of heat, so the temperature rises rapidly, which is the first stage. When the fluid temperature near the wall rises to the pseudo-critical temperature, cp increases rapidly. The working fluid then absorbs more heat, so the wall temperature does not increase as quickly and the fluid temperature continues to increase. In this stage, fluid temperature distribution is more uniform in the radial direction and close to the pseudo-critical temperature as shown in Fig. 10, so the average fluid cp is higher. As the overall fluid temperature in the tube increases, the proportion of fluid in the high cp region changes with the average cp first increasing and then decreasing. This is the second stage. When the centerline fluid temperature exceeds the pseudo-critical temperature, cp and the fluid thermal conductivity rapidly decrease, which reduces the heat transfer coefficient in the third stage.

Fig. 10
Radial fluid temperature distributions along the Y axis As mentioned earlier, the flow in smaller inner diameter tubes has higher wall and fluid temperatures. In the 0.5 mm inner diameter serpentine tube, the fluid temperature near the exit exceeds the pseudocritical temperature, so this smallest diameter tube has all three temperature distribution stages.
Therefore, the heat flux should be matched to the inner diameter and length of the serpentine tube when designing microchannel heat exchangers.

Effect of the geometric parameters and flow direction on the average convective heat transfer coefficient
Constant mass flow rate The effects of the geometric parameters and flow direction on the average convective heat transfer coefficient are shown in Fig. 11. The average convective heat transfer coefficient decreases significantly with increasing inner diameter when the tube pitch is held constant. The increasing inner diameter will reduces the velocity with a constant mass flow rate, and the heat transfer coefficient decreases with the flow velocity. At the same time, the average convective heat transfer coefficient decreases with increasing tube pitch when the inner diameter is held constant because a larger tube pitch reduces the centrifugal force which then reduces the turbulence intensity and the heat transfer coefficient.

Fig. 11 Heat transfer coefficients for different flow directions
In addition, the results show that the flow direction only affects the average convective heat transfer coefficient for inner diameters larger than 0.5 mm. For inner diameters larger than 0.5 mm, the average convective heat transfer coefficient is higher for downward ScCO2 flow. This shows that buoyancy has a major effect on the heat transfer coefficient for ScCO2 flow in larger diameter tubes. The CO2 density decreases sharply near the pseudo-critical point, which leads to buoyancy effects and flow acceleration. With upward flow, the buoyancy and flow acceleration are in the same direction and reduce the velocity gradient and shear stress between the boundary layer and main stream which reduces the heat transfer coefficient. With downward flow, the buoyancy reduces the influence of the flow acceleration, increases the turbulent shear stress and the velocity gradient near the wall, increases the turbulence intensity, and enhances the heat transfer coefficient. As shown in Fig. 12, the wall temperature is lower for downward ScCO2 flows than for upward flows. Bo increases with increasing inner diameter for the same inlet mass flow rate as shown in Fig. 13 with the buoyancy effect known to be important for larger Bo.  shows an enlarged graph of the curves for d=0.5 mm shown in Fig. 13 (a). The results show that the flow direction has no effect on the heat transfer coefficient for the 0.5 mm inner diameter. Thus, the buoyancy effect can be ignored for Bo less than 1×10 -7 in serpentine tubes. For the inner diameters of 1 mm, 1.5 mm and 2 mm, Bo＞3.6×10 -6 and the buoyancy effect cannot be ignored. In these cases, the turbulent shear stress and velocity gradient near the wall for the downward flow are greater than for the upward flow because of the buoyancy effect and the turbulence intensity for the downward flow is greater than for the upward flow. Therefore, the heat transfer coefficient for downward flow is larger than for upward flow as shown in Fig. 11. The Bo distribution has a peak because the CO2 density changes dramatically with small variations in the CO2 temperature near the pseudo-critical temperature. These large radial density variations then lead to large buoyancy effects. Fig. 14, the inner diameter, tube pitch and flow direction all significantly affect the average convective heat transfer coefficient. The average convective heat transfer coefficient in the serpentine tube first increases and then decreases with increasing inner diameter. Thus, the inner diameter has an optimal value which gives the maximum average convective heat transfer coefficient. Moreover, this optimal heat transfer coefficient occurred at the same diameter for the three tube pitches used in this study. This phenomenon is related to the relative magnitudes of the influence of the heat exchange area and the radial temperature difference. A smaller radius leads to a smaller radial temperature difference, a more uniform fluid temperature distribution across the cross-section, and a higher proportion of the working fluid in the high specific heat region which increases the heat transfer coefficient. However, a small inner diameter also reduces the heat transfer area, which reduces the heat transfer rates. For small inner diameters, the heat transfer area is most important in controlling the heat transfer, so the heat transfer increases with increasing inner diameter. As the diameter increases, the heat transfer area changes have less effect while the radial temperature difference has a stronger effect and the heat transfer decreases with increasing diameter. Therefore, the average convective heat transfer coefficient in serpentine tubes has a maximum the increasing inner diameter. For a constant inlet Reynolds number, the average convective heat transfer coefficient in the serpentine tube decreases with increasing tube pitch with a higher average convective heat transfer coefficient for downward ScCO2 flow.

CONCLUSIONS
The flow and heat transfer for heating of ScCO2 were modeled in various serpentine tubes in this study. Twelve three-dimensional geometries were generated with various inner diameters (0.5, 1.0, 1.5, and 2.0 mm) and tube pitches (4, 8 and 16 mm). The models were used to analyze the effects of the geometry and flow direction for constant inlet mass flow rate and constant inlet Reynolds number boundary conditions.
No heat transfer deterioration was observed for upward flow in vertical serpentine tubes in this study due to the effects of the centrifugal force on the flow. When the mass flow rate was held constant, the outer tube wall temperature increases with increasing tube pitch for a constant inner diameter due to the reduction of the centrifugal force and the turbulence intensity. A smaller inner diameter and smaller tube pitch make the radial temperature distributions more uniform. The average convective heat transfer coefficient decreases as both the inner diameter and the tube pitch increase. The heat transfer coefficients were higher for downward flow than for upward flow due to the effect of buoyancy on the turbulent intensity for inner diameters larger than 0.5 mm. When the inlet Reynolds number was held constant, the outer wall temperature and the centerline fluid temperature both increase with decreasing inner diameter because of the reduced heat transfer area. The average convective heat transfer coefficient first increases and then decreases with increasing inner diameter due to the relative magnitudes of the influence of the heat transfer area and the radial temperature difference. The average convective heat transfer coefficient decreases with increasing tube pitch. In addition, the heat transfer coefficients were higher with downward ScCO2 flow than with upward flow, and the buoyancy effect can be ignored for Bo＜1×10-7.
The present work will facilitate optimization of serpentine tube diameters, pitches and lengths for heat exchangers to improve heat exchange efficiencies and system performance.