Abstract

Supercritical carbon dioxide (sCO2) can be utilized as a working fluid in various thermal systems including large-scale power cycles; portable power production units, centralized coolant systems, and standalone cooling devices. However, the lack of accurate prediction tools such as heat transfer coefficient correlations and insufficient research studies about the mechanisms controlling heat transfer processes are hindering its practical realization for key energy and cooling systems. The overall objective of this study is to extend fundamental knowledge about heat transfer and fluid dynamic processes in conduits pertinent to sCO2 with an emphasis on flow direction and inclination effects. This paper presents the study on effects of gravity, buoyancy on sCO2 flow at temperatures near and away from the pseudo-critical temperature. The experimental setup consists of a high temperature and pressure sCO2 heat transfer loop and flow testing facility. Recently, researched sCO2 heat exchangers can have tubes oriented at different angles such as 45 deg or 90 deg to horizontal. For optimized design of efficient and cost-effective turbomachinery components utilizing sCO2 as the heat transfer fluid, an understanding of convective heat transfer inside a tube/pipe is equally as important as external heat transfer. This paper presents an experimental and numerical study on sCO2 heat transfer at various inclinations with angles ranging from 0 deg (horizontal) to 90 deg (vertical) along with upward and downward flow direction with different inlet temperatures. Thermocouple-based temperature measurement is utilized at multiple locations within the test section axially and circumferentially to study the temperature distributions on the tube surface. Computational fluid dynamics (CFD) simulations have been performed using ANSYS Fluent to complement experimental data. The CFD and experiment have been analyzed against the well-known Gnielinski Nusselt number correlation.

1 Introduction

Supercritical CO2 heat transfer has gained a lot of recent attention due to its wide range of potential applications in areas, including, concentrated solar power plants [1], refrigeration [2,3], waste heat recovery [4,5], thermonuclear power generation [6], and air-conditioning systems [7,8]. The critical pressure and temperature for CO2 are 7.38 MPa and 304 K, respectively. Near the critical point, the fluid undergoes rapid changes in its thermophysical properties which makes it important to understand its heat transfer behavior if it is to be used as a working fluid in thermal engineering applications. There have been numerous experimental and numerical studies pertaining to sCO2 heat transfer in circular cross section tubes. Most of the literature is directed toward horizontal or vertical flow direction and there is a lack of knowledge when it comes to studying heat transfer for flow at an arbitrary inclination angle. There has been some numerical work for the latter [8], but no experimental work has been found within the public domain.

As multiple authors have shown, there can be significant influence of property variations on sCO2 heat transfer close to the pseudo-critical temperature [915]. Near the critical point of CO2, Schnurr [9] reported circumferential variation in temperature at a cross section in a horizontal tube for sCO2. Adebiyi and Hall [10] explained that occurrence of natural convective current from bottom surface to top surface of the horizontal tube is the reason behind enhanced heat transfer at the bottom surface. This phenomenon is represented by a sketch as shown in Fig. 1. This causes higher temperature at the top surface compared to the bottom surface of the tube. Liao and Zhao [16] reported that this effect of buoyancy on heat transfer is directly related to size of the tubes. This is mainly due to the Richardson number, which is a measure of natural convection and a ratio of Grashof number to Reynolds number squared.

Fig. 1
Schematic showing formation of natural convection currents in a cross section of a heated circular tube at inclination angles not equal to 90 deg
Fig. 1
Schematic showing formation of natural convection currents in a cross section of a heated circular tube at inclination angles not equal to 90 deg
Close modal

In the vertical orientation, sCO2 heat transfer in tubes is highly dependent on flow direction. The effects of buoyancy on vertical flow heat transfer under turbulent conditions are explained by Jackson and Hall [17]. Heat transfer at the wall is directly related to variation in generation of wall shear stress. For upward direction flow, buoyancy forces are opposing wall shear stress, decreasing turbulence and hence, result in decreasing heat transfer. On the other hand, for downward direction flow, buoyancy forces increase wall shear stress, generating increased turbulence and hence, increasing heat transfer at the wall.

This paper presents experimental and numerical analysis on sCO2 heat transfer in a heated tube at five different inclinations; 0 deg, 30 deg, 45 deg, 60 deg, and 90 deg. Both upward and downward direction flow cases are performed here. Two inlet temperatures are considered here; one near the pseudo-critical temperature, and one away, to study their effects on heat transfer.

2 Experimental Methodology

2.1 Description of Experimental Rig.

A schematic of the experimental rig can be seen in Fig. 2. The test loop is setup in a closed-loop configuration to conserve usage of CO2 and is designed for maximum operating pressure and temperature of 100 bar and 150 °C. As seen from the right image in Fig. 3, the test section is located within a structural frame to facilitate variation of inclination by adjusting winches, hooked to the test section frame. Before filling the loop with CO2, the loop is vacuumed to remove any air, ensuring close to pure concentrations of CO2 for testing. The loop is also equipped with an exhaust system to safely discharge CO2 from the loop. A CO2 cylinder and booster pump are used in series to pressurize the loop to the desired pressure. After the filling process, the booster pump is disconnected from the loop using an ON/OFF valve. A micropump gear pump is used to circulate sCO2 through the loop, with a maximum operating pressure of 105 bar. A buffer tank is used downstream of the gear pump to absorb any possible pressure fluctuations and deliver steady mass flowrate to the test section. A Coriolis mass flowrate meter suitable for CO2 flows in gas, liquid, and critical phases is used to measure mass flowrate. Inlet temperature to the test section is controlled using a preheater made from wrapped-around electric rope heaters and voltage controlling variacs supplying constant voltage.

Fig. 2
Schematic of closed-loop experimental setup to study sCO2 heat transfer
Fig. 2
Schematic of closed-loop experimental setup to study sCO2 heat transfer
Close modal
Fig. 3
(a) Photograph of the test section inclined at 45 deg , (b) 3D cad model of the frame holding test section, and (c) winches
Fig. 3
(a) Photograph of the test section inclined at 45 deg , (b) 3D cad model of the frame holding test section, and (c) winches
Close modal

Electric heating of the test section is achieved using a DC power supply, which is connected to the test section tube using machined copper busbars. These machined copper busbars ensure low thermal contact resistance achieved from good contact with test section tube. A busbar, which is in loose contact with the tube, might create local spikes in resistance, resulting in heat generation, making it difficult to account for heat transferred to the flow. Seamless stainless steel tubing of outer diameter 12.7 mm is used in the test section and tubing connections throughout the experimental loop. The test section is made of ASTM A213 SS316 tube with outer diameter of 12.7 mm and inner diameter of 10.92 mm, enclosed in a polyvinyl chloride pipe and insulated using mineral wool insulation. A water cooler is used downstream of the test section before the recirculating pump to prevent from exceeding temperature above working limit of the pump, as well as to achieve overall steady-state during testing. The water cooler consists of a tube array immersed in a cold water drum where the water temperature is controlled and recycled using a chiller.

Inlet flow temperature is measured using two different instruments to add redundancy in measurement. These consist of a K-type thermocouple probe and a RTD probe inserted at the same axial location using a cross-union fitting. The same setup is used to measure outlet flow temperature as well. Static pressure is measured at the inlet using a pressure transducer. The test section includes an unheated section with L/Dint = 74, heated section with L/Dint = 123 and unheated section L/Dint = 74 before measuring outlet flow temperature and the busbar sections L/Dint = 1.15. These values of L/Dint ratios are chosen based on conventional values necessary to fully develop the hydrodynamic boundary layer and thermal boundary layer. The data reduction process takes into account outside wall temperatures in the test section, which are measured using T-type thermocouples instrumented on the tube surface. At each axial location, four surface thermocouples are attached to the surface at 90 deg tangential interval to measure temperatures at top, bottom, right, and left wall locations (as shown in right side of Figs. 2 and 4(b)). Electric heating is accounted by measuring current provided by the DC power supply and voltage measured between two busbars.

Fig. 4
(a) Test section and (b) a domain of heat transfer coefficient calculation
Fig. 4
(a) Test section and (b) a domain of heat transfer coefficient calculation
Close modal

2.2 Heat Loss Tests.

Even though the stainless steel tube is insulated using mineral wool insulation, some heat loss through the insulation to the surrounding will occur. Hence, heat loss calculations are important to calculate heat transfer coefficients correctly. Heat loss from the tube to ambient is a function of temperature difference between tube external wall and ambient (ΔTw−amb), outer tube area and area-specific equivalent conductance (Uloss (W/m2K)) between outer wall of tube and ambient as shown in Eq. (1). The goal of heat loss tests was to calculate this equivalent conductance, Uloss, since the other two parameters are known. ΔTw−amb is measured using thermocouples instrumented on outer wall of tube and ambient. To calculate Uloss, the test loop is vacuumed so that the inside wall of the tube does not experience any surface heat convection. Here, it is assumed that the axial conduction from the heated section of the tube to unheated section of the tube is negligible. That is why Eq. (1) holds true when the test loop is vacuumed
ΔQelectric,local=ΔQloss,radial,local=UlossΔTwambLlocalπDext
(1)

Three heat loss tests are performed with three different values of electric power to cover the range of ΔTw −amb expected from flow experiments. Uloss is calculated for each case and plotted against ΔTw −amb as show in Fig. 5 to obtain an expression for Uloss. Heat loss is negligible, with a maximum heat loss percentage of 0.056% calculated between the electrical power and the heat added to the test section, for all experimental cases. Heat loss tests to calculate the equivalence conductance were performed for completeness. This plot is shown including the expression for Uloss as a function of ΔTw −amb in Fig. 5.

Fig. 5
Heat loss coefficient obtained from no-flow experimental tests
Fig. 5
Heat loss coefficient obtained from no-flow experimental tests
Close modal

2.3 Data Reduction.

The test section is divided into seven zones as shown in Fig. 4(a) to account for an accurate energy balance. These zones include energy balance based on difference in bulk enthalpy between outlet and inlet. Figure 4(b) also shows all possible pathways of heat transfer in a domain, where Nusselt number is calculated. The approach to calculate local bulk enthalpy includes axial energy balance marching considering heat transferred to CO2. This energy marching starts at the inlet where inlet enthalpy is calculated using REFPROP [18] from measured bulk temperature and pressure. Figure 4(b) shows the domain for calculation of heat transfer coefficient along with pathways of heat transfer to CO2. Each heat transfer domain consists of 8.7 × Dint length around a thermocouple station and 90 deg sector (quadrant) of annular cross section of tube, each for four circumferential locations considered—top, bottom, left, and right. Data reduction process to calculate local bulk temperature and local Nusselt numbers follows the same equations as reported in the previous study [19].

Local electrical heat addition necessary for heat balance is proportional to the local length of the domain of heat transfer coefficient calculation and total length of heated tubing as shown in the following equation:
ΔQelectric,local=Qelectric,total×ΔLlocalLtotal
(2)
The calculated volumetric heat generation is used to estimate internal wall temperature, along with the thermal conductivity of stainless steel. Volumetric heat generation (qvl) is calculated using total electrical heat generation (Qelectric,total) and solid metal volume of the heated portion of the tube, as shown in Eq. (3). Thermal conductivity of stainless steel being a function of temperature is calculated using Eq. (4) [20]. Internal wall temperature is then obtained using Eq. (5).
qvl=POWERLtotal×Acs
(3)
kSS=14.6+1.27×102×Text
(4)
Tw,int=Tw,ext+qvl/(4kSS)[rext2rint2]qvl/(2kSS)(rext)2ln(rext/rint)
(5)
Heat loss to the ambient in the radial direction is calculated using equivalent loss coefficient Uloss and temperature difference between external wall and ambient (ΔTw −amb), as shown in Eq. (6). Area of heat loss considered is the local external surface area of the tube (Llocal× π × Dext)
ΔQloss,radial=Uloss×ΔTwamb×ΔLlocal×π×Dext
(6)

Figure 4(b) shows the domain for calculation of heat transfer coefficient along with pathways of heat transfer to CO2. Each heat transfer domain consists of a length 8.7 × Dint around a thermocouple station and a 90 deg sector (quadrant) of annular cross section of tube, each for four circumferential locations considered (top, bottom, left, and right.)

Heat is convected to CO2; hence, electrical heat added minus heat lost to surrounding, as shown in the following equation:
ΔQCO2,local=ΔQelectric,localΔQloss,radial
(7)
Bulk flow enthalpy (hb,i) at a specific axial station (i) is considered common for all quadrants. hb,i is calculated by adding Δhb,i to the enthalpy calculated at the previous axial station (hb,i −1), as shown by Eq. (8). Δhb,i is calculated using Eq. (9)
hb,i=hb,i1+Δhb,i
(8)
where
Δhb,i=Quadrant1Quadrant4ΔQCO2,localm˙
(9)
Bulk flow temperature Tbulk,i at a specific axial station (i) is calculated using REFPROP [18] with local enthalpy and pressure as inputs
Tb,i=fn(hb,i,P)
(10)
Temperature difference between the bulk flow and the wall is used to calculate heat transfer coefficient as follows:
HTClocal=ΔQCO2,local/(m˙×ΔTwb,local×Alocal)
(11)
CO2 bulk flow thermal conductivity, kb, is calculated using REFPROP with bulk flow temperature and pressure as inputs. This value of CO2 thermal conductivity is used in calculation of local Nusselt number as shown in the following equation:
Nulocal=HTClocal×Dint/kb
(12)
Gnielinski correlation Nusselt number [21] (Eq. (13)) is used for comparison with experimentally calculated Nusselt number
NuGn=(f/8)(Reb1000)Prb1+12.7(f/8)0.5(Prb2/31)
(13)
Grashof number, Reynolds number, and Richardson numbers are defined in the following equations, respectively:
Gr=gρb(ρbρw)Dint3μb2
(14)
Reb=4m˙πμbDint
(15)
Ri=GrReb2=gρb(ρbρw)Dint516m˙2
(16)

Since mass flowrate for all the cases presented here is constant, Reynolds number equation with mass flowrate is adopted here.

2.4 Uncertainty Analysis.

The uncertainty analysis is based on methods described by Kline [22] and the Test Uncertainty Standard PTC 19.1 by the American Society of Mechanical Engineers (ASME) [23]. Uncertainties of different measured parameters are listed in Table 1. Table 2 shows a summary of uncertainty calculated for Nusselt numbers with a maximum uncertainty of 6%.

Table 1

Uncertainties of measured parameters

ParameterInstrument uncertainty
Test section power±1.63%
Inlet pressure±1.9%
Temperatures1.1 °C or 0.75% whichever is larger
Mass flow rate±0.98%
ParameterInstrument uncertainty
Test section power±1.63%
Inlet pressure±1.9%
Temperatures1.1 °C or 0.75% whichever is larger
Mass flow rate±0.98%
Table 2

Nusselt number uncertainty summary

Maximum uncertainty6%
Minimum uncertainty2.5%
Average uncertainty4%
Maximum uncertainty6%
Minimum uncertainty2.5%
Average uncertainty4%

2.5 Air Validation.

Due to the unconventional nature of the experimental setup described here, it is important to validate the correctness of experimental results and the data reduction process obtained from this setup. For this paper, validation experiments with high pressure air were conducted to establish the baseline confidence interval for the sCO2 tests and the data reduction procedure.

Figure 6(a) shows comparison of experimentally obtained Nusselt numbers and Gnielinski Nusselt number for a case with the mentioned testing conditions. The match for this case was within ±3%. Overall, for all air tests, the Nusselt number calculated from experiments fit well with predicted Nusselt numbers from Gnielinski correlations. This is plotted in Fig. 6. The largest error between the predicted Nusselt numbers and calculated values was less than ±6% for the range of Reynolds number tested. This is within maximum uncertainty observed in the experimentally calculated Nusselt numbers.

Fig. 6
Heat transfer validation results using air: (a) Nusselt numbers plotted at axial locations for a case and (b) Nusselt number from all cases plotted against Reynolds number (Error bars and dotted lines are ±6% of Gnielinski Nusselt number)
Fig. 6
Heat transfer validation results using air: (a) Nusselt numbers plotted at axial locations for a case and (b) Nusselt number from all cases plotted against Reynolds number (Error bars and dotted lines are ±6% of Gnielinski Nusselt number)
Close modal

3 Numerical Methodology

3.1 Geometry and Numerical Modeling.

The companion numerical simulations performed in this study aim to model the experimental configuration as closely as possible; this includes geometry and boundary conditions. To allow for direct comparison to experimental results, surface temperature probes and bulk temperatures are monitored at identical streamwise, circumferential locations as those of the physical experiment and the same volumetric heat generation for the heated section (zone 4) of the tube as shown in Fig. 4. Solid domain pipe dimensions are identical to that of the experiment and are given in Sec. 2. Numerical boundary conditions for the companion simulations used a pressure inlet (based on case experimental values of pressure and bulk temperature), mass flow exit, volumetric heat generation for the solid domain, and an external wall temperature-dependent convective boundary condition for the solid domain, as shown in Table 3 This temperature-dependent convective boundary condition uses the experimental heat loss relation from Fig. 5, to model experimental heat loss.

Table 3

Experimental and CFD conditions

Inlet pressure (MPa)8.0
Mass flow rate (Kg/s)8.33 × 103
Inlet temperature (pseudo-critical) (°C)35
Inlet temperature (°C)55
Volumetric heat generation (W/m3)1.45 × 107
Angles of inclination (°)0 (horizontal), 30, 40, 60, 90 (vertical)
Flow directionsUpward and downward
Inlet pressure (MPa)8.0
Mass flow rate (Kg/s)8.33 × 103
Inlet temperature (pseudo-critical) (°C)35
Inlet temperature (°C)55
Volumetric heat generation (W/m3)1.45 × 107
Angles of inclination (°)0 (horizontal), 30, 40, 60, 90 (vertical)
Flow directionsUpward and downward

The simulations for this study were run with the finite volume commercial solver ANSYS Fluent. Steady-state conjugate heat transfer simulations for sCO2 were run with double precision and 2nd order upwind discretization schemes. Fluent's pressure-based solver was used with the fully coupled scheme for pressure–velocity coupling, with the PRESTO! scheme used for pressure interpolation. With strong expected buoyancy effects, the PRESTO! scheme was chosen for its ability to model such flows. The current problem is treated as a single-phase flow and thus can be solved with the general fluid and heat transfer equations. The general equations are given as:

Mass conservation equation
(ρui)xi=0
(17)
Momentum conservation equation
(ρuiuj)xj=Pxi+xj[(μ+μt)(uixj+ujxi)23(μ+μt)ukxkδi,j]+ρgi
Energy conservation equation
(ρTuj)xj=xj[(kCp+μtPrt)Txj]
(18)
Heat conduction equation for solid domain
·(kssT)=0
(19)

For this study, the k − ω Shear Stress Transport turbulence model [24] was used based on its documented success in challenging, buoyancy-impacted sCO2 flows [2527]. Since this study focuses on inlet bulk temperatures close to the critical point where thermophysical property variations are considerably high, CO2 is modeled as a real gas with variable properties. Property lookup tables are generated within Fluent's UDF library using the NIST Database Version 9.1 (REFPROP v9.1). The temperature resolution of the generated lookup table is 0.065 °C. Solid domain material properties for 316 Stainless Steel used a constant density and temperature dependent functions for specific heat and thermal conductivity.

Numerical meshes for the fluid and solid domains were generated within ANSYS workbench and consisted of unstructured tetrahedral and hexahedral cells, respectively. The first cell height in the fluid domain prism mesh was set to maintain wall y+ < 1 throughout the entire domain. Twenty-two prism layers were used in the inflation layer. Due to the presence of buoyancy-driven secondary flows within the heated section of the domain, additional volumetric refinement was performed (zones 2 through 6) to improve resolution in the core flow area. The solid domain mesh maintained similar cell base size to that of the fluid domain for consistency and robustness.

The pipe domain mesh maintained 10 cells in the radial direction, with smaller cell heights biased toward the inner and outer wall boundaries.

Mesh independence was verified by the evaluation of external wall temperature distributions, across four numerical meshes of varying count (7, 9, 13.7, 16, and 20 million). Figure 7(a) shows the minimum, maximum, and average temperature difference of each mesh compared to the finest mesh of 16 million. “Minimum” and “Maximum” refer to the highest and lowest absolute temperature deviations between the respective meshes. Mesh independence is reached with a cell count of 13.7 million, where the deviation between it and the 20 million mesh is less than 1 °C. Note that all meshes evaluated retained a y+< 1. Coarsening/refinement was accomplished via base size adjustment and volumetric refinement within the heated portion of the domains. All simulation results presented herein use the 20 million count mesh.

Fig. 7
Mesh convergence study
Fig. 7
Mesh convergence study
Close modal

The CFD data are validated using temperature, by comparing the zone 4 temperatures for the 35 °C and 55 °C to the experimental external wall temperatures. In order to obtain the wall temperature, the simulation is ran until the wall temperatures reached steady-state and convergence is achieved. The results of the simulation study are presented in Figs. 8 and 9 along with the experimental data, the CFD effectively predicts the wall temperature trend although for the case of 35 °C in Fig. 8, the trend matches for both CFD and experiment in the axial direction, but there is mismatch in the radial direction, particularly at the upper wall of the pipe. With the 35 °C case being near the pseudo-critical temperature, the heat transfer prediction becomes challenging, as is noted in other numerical works in the literature. The 55 °C inlet temperature case shows excellent match with the experimental results.

Fig. 8
Zone 4 wall temperature validation of CFD with experiment for 35 °C horizontal flow case
Fig. 8
Zone 4 wall temperature validation of CFD with experiment for 35 °C horizontal flow case
Close modal
Fig. 9
Zone 4 wall temperature validation of CFD with experiment for 55 °C horizontal flow case
Fig. 9
Zone 4 wall temperature validation of CFD with experiment for 55 °C horizontal flow case
Close modal

4 Results and Discussion

This study analyzed the effect of inclination on the heat transfer properties of sCO2 flowing in a tube at an operating pressure of 8.0 MPa, for a mass flux of 90 kg/m2s and heat flux of 18 kW/m2 and 12 kW/m2, with inlet bulk temperatures of 35 °C and 55 °C, respectively. The effect of parameters, such as varying pressure and mass flux, has been studied by Pidaparti et al., which showed that when the bulk temperature is slightly lower than the pseudo-critical temperature, the heat transfer coefficient is highest for lower operating pressures and when the pseudo-critical temperature is higher than the bulk temperature, the heat transfer coefficient is higher for higher operating pressures [13]. The results will show the effect of inclination on the heat transfer properties of the flow. The plots and images shown below are for parameters in zone 4, from Fig. 4, zone 4 consists of 8 thermocouple stations numbered 5–12. The plots are for stations 6–11 while the contour profiles are for stations 5–11.

4.1 Effect of Inlet Bulk and Wall Temperature.

The results of the 35 °C and 55 °C cases for horizontal and 90 deg flow direction show the effect of inlet bulk temperature on the flow heat transfer correlation. At 35 °C, which is approximately the pseudocritical temperature at 8 MPa, for the horizontal case, the wall temperature at the top is significantly higher than the wall temperature at the side and the bottom, with the bottom showing the lowest temperature in the radial direction of the tube. This trend is consistent across all stations as shown in Fig. 10.

Fig. 10
Experimental external wall temperatures and bulk temperatures for horizontal flow cases: Tin = 35 °C and Tin = 55 °C
Fig. 10
Experimental external wall temperatures and bulk temperatures for horizontal flow cases: Tin = 35 °C and Tin = 55 °C
Close modal

This indicates that there is temperature variation of wall temperature in the tube for the horizontal case, which becomes less prevalent as the tube is inclined. This circumferential variation of wall temperature for the horizontal case is due to the density difference in the tube. The density variation in the tube creates a region where the high density fluid settles at the bottom of the fluid domain, as indicated in Fig. 1. The phenomenon causes buoyancy and flow acceleration due to the difference in density of the bulk flow and the near wall flow, hence creating turbulent shear stress at the higher density region (bottom of the tube), which locally enhances the fluid's thermal conductivity and heat transfer; this is verified in Figs. 11 and 12, which shows that the bottom of the tube, especially for the horizontal case, has the highest Nusselt number. The phenomenal effect of buoyancy and gravity defining the thermophysical properties of the flow is further discussed in the next subsection, effect of flow direction and buoyancy.

Fig. 11
Experimental Nusselt number plotted against axial locations (z/Dint) in heated section for 35 °C inlet temperature at different inclinations
Fig. 11
Experimental Nusselt number plotted against axial locations (z/Dint) in heated section for 35 °C inlet temperature at different inclinations
Close modal
Fig. 12
Experimental Nusselt number plotted against axial locations (z/Dint) in heated section for 55 °C inlet temperature at different inclinations
Fig. 12
Experimental Nusselt number plotted against axial locations (z/Dint) in heated section for 55 °C inlet temperature at different inclinations
Close modal

The inlet bulk temperature affects the heat transfer properties for the 35 °C horizontal case; the experimental Nusselt number increases as the bulk temperature increases, which is opposite of the Gnielinski correlation Nusselt number as shown in Fig. 11 for the horizontal flow case. The Nusselt number trend for inlet bulk temperature 55 °C for the horizontal case showed better matching and prediction for the sidewall experimental Nusselt number and trend in the same direction. This further shows that the Gnielinski correlation can be used in cautiously predicting the Nusselt number for inlet temperatures which are significantly higher than the pseudo-critical temperatures. Figure 10 shows the difference between the top and bottom wall temperatures and can be seen to decrease as the inlet bulk temperature increases. Due to the difference in the bulk density and the near wall density, decreasing buoyancy effects are further observed compared to the 35 °C inlet bulk temperature case.

4.2 Effect of Flow Direction and Buoyancy.

Figure 10 shows the temperature profile for the horizontal cases of 35 °C and 55 °C, there is notable difference between the wall and bulk temperatures, which can be explained as the effect of radial convection as heat is transferred from the tube to the fluid, as the inlet temperature is increased, the effect of buoyancy is diminished. The plots show the external wall and bulk temperatures at locations inside the heated section. There is a temperature gradient at each z/Dint location, which is due to radial density variation. There is also an increase in the bulk and wall temperatures with the top wall temperature being the highest. As the inlet bulk temperature is increased from 35 °C to 55 °C, there is significant increase in the wall and bulk temperatures in the axial direction, although the temperature gradient at each z/Dint location (radial direction) does not reflect the increase; from Fig. 10, the temperatures at the wall and bulk maintain about the same trendline in the axial direction; hence, the difference between the wall temperatures and the difference in any wall temperature and at each z/Dint location is approximately the same. The temperature difference between the top and bottom wall is decreasing as from ΔT = 22 °C for the 35 °C horizontal case, to ΔT = 16 °C for the 55 °C horizontal case; this shows that with an increase the inlet bulk temperature away from the pseudo-critical temperature, the radial temperature gradient of the wall decreases while the temperature difference between the sidewall and the bulk temperature at the z/Dint location increases.

The Nusselt number plots in Figs. 11 and 12 show the trend of Nusselt number calculated from the internal wall temperature using heat transfer coefficient and Gnielinski correlation, the Gnielinski Nusselt number follows a downward trend for both cases, which indicates that as the bulk temperature is increasing in the axial direction, the Nusselt number decreases, which is relatively similar to the experimental Nusselt number trend for the 55 °C case. The Nusselt number plot in Fig. 11 for the 35 °C case shows the experimental Nusselt number increasing in the axial direction while the Gnielinski Nusselt number decreases, however.

This singularity of trend difference in the plot is as a result of the different inlet bulk temperature and its thermophysical properties. In Fig. 10, it is shown that as inlet bulk temperature increases to be sufficiently higher than the pseudo-critical temperature, there is an increase in ΔT between the wall and bulk temperatures. Analyzing Fig. 13, we see the effect of inlet bulk temperature on thermal conductivity, which is inversely proportional to the experimental Nusselt number. The rate of decrease of the thermal conductivity for inlet bulk temperature of 35 °C closer to the pseudo-critical temperature is very rapid when compared to that of the 55 °C case; it is interesting to note that this, coupled with the high thermal capacitance of the fluid at these near pseudo-critical bulk temperatures (less than 40 °C), greatly affects the upward and downward trend of the experimental Nusselt number in Figs. 11 and 12.

Fig. 13
Axial variation in thermal conductivity of CO2 for horizontal flow cases
Fig. 13
Axial variation in thermal conductivity of CO2 for horizontal flow cases
Close modal

The effect of buoyancy, which is largely due to density effects, can be reduced by inclination, which affects the turbulent shear stress at the wall of the tube. The effect of buoyancy can be characterized by a nondimensional parameter Grashof and/or Richardson number, for flow in a pipe or tube. In order to investigate the effect of buoyancy on the flow, we analyzed the Nusselt number plot and CFD contours for horizontal, 45 deg and 90 deg case for inlet bulk temperatures of 35 °C and 55 °C. Buoyancy forces act upwards and opposite to gravitational forces. Figure 14 shows the Richardson number obtained from experiment, plotted against the axial distance. For the case of 45 deg downward flow, the effect of buoyancy, which is characterized by Richardson number, is decreasing rapidly in the axial direction for the 35 deg case. The 55 deg case shows a decrease in the Richardson number in the axial direction, but it is evident that for inlet bulk temperatures closer to the pseudo-critical temperature, radial density variation at each z/Dint is predominant in the flow.

Fig. 14
Richardson numbers variation calculated from experimental data plotted against axial locations (z/Dint) in heated section for 35 °C and 55 °C inlet temperature at different inclinations and flow directions
Fig. 14
Richardson numbers variation calculated from experimental data plotted against axial locations (z/Dint) in heated section for 35 °C and 55 °C inlet temperature at different inclinations and flow directions
Close modal

In order to visualize the effect of buoyancy in the flow, CFD was performed. It is observed in Figs. 15 and 16 that there is density variation at each z/Dint location with the heavier fluid at the bottom of the tube. The temperature at the tube inlet is uniform, but as the flow enters the heated section and heat is added, the fluid develops a radial temperature gradient with peak temperature at the top wall location. This creates the radial density variation shown in the contour plot, and hence the buoyancy forces deteriorate as the fluid local bulk temperature is increased in the axial direction and as the flow direction changes. Buoyant forces act on the fluid creating an “m”-shaped density field with the secondary flow fluid at the top of the tube, and the heavier fluid at the bottom. Figures 11 and 12 also validate the density contour plots; the effect of buoyancy can be seen on the separation between the Nusselt number at the top and bottom wall. As the flow direction changes, the separation between the top and bottom reduces until the flow experiences negligible buoyant forces.

Fig. 15
CFD Contour plot of density and velocity profile for upward flow at Tin = 35 °C showing cross section cuts at axial locations (z/Dint) in heated section. Density contour: (a) 90 deg (b) Horizontal flow (c) 45 deg. Velocity contour: (e) Horizontal flow (f) 45 deg. Flow direction and axial location (d).
Fig. 15
CFD Contour plot of density and velocity profile for upward flow at Tin = 35 °C showing cross section cuts at axial locations (z/Dint) in heated section. Density contour: (a) 90 deg (b) Horizontal flow (c) 45 deg. Velocity contour: (e) Horizontal flow (f) 45 deg. Flow direction and axial location (d).
Close modal
Fig. 16
CFD Contour plot of density and velocity profile for downward flow at Tin = 35 °C showing cross section cuts at axial locations (z/Dint) in heated section. Density contour: (a) 90 deg (b) 45 deg. Velocity contour: (d) 45 deg. Flow direction and axial location (c).
Fig. 16
CFD Contour plot of density and velocity profile for downward flow at Tin = 35 °C showing cross section cuts at axial locations (z/Dint) in heated section. Density contour: (a) 90 deg (b) 45 deg. Velocity contour: (d) 45 deg. Flow direction and axial location (c).
Close modal

In Fig. 16(d), contour fields of in-plane velocity magnitude show slightly increased boundary layer thickness when compared to velocity contours of horizontal flow in Fig. 15(e). One can also note lower near-wall velocity magnitudes when compared to horizontal, or 45 deg upward flow. As the angle of inclination is increased for the case of downward flow, the turbulent shear stress is enhanced by the buoyant forces acting on the fluid which then creates localized drops in the wall temperature, thereby increasing the heat transfer and the Nusselt number. In the case of upward flow, as the angle of inclination is increased, buoyant forces reduce to become relatively negligible at the limiting case of 90 deg. Both upward and downward 90 deg cases are highlighted by the lack of large scale secondary buoyancy driven structures. Velocity contours are not shown for these cases, as the in-plane velocity magnitudes are an order of magnitude smaller, than that of the horizontal and 45 deg inclined cases. The density fields seen for Figs. 15(a) and 16(a) are axisymmetric, which support the observed experimental station wall temperature uniformities.

5 Conclusions

This work was performed to understand the fluid dynamic and heat transfer properties associated with high pressure and varying inlet temperature flow of sCO2 in a tube. Experimental and CFD studies were performed at a pressure of 80 MPa, mass flux of 90 Kg/m2s with inlet bulk temperatures of 35 °C and 55 °C. The study showed that the flow direction and inlet temperature are major parameters affecting the heat transfer properties of sCO2. Nondimensional parameters such as Nusselt and Richardson number are plotted from experimental data to show the effect of the varying parameters on heat transfer and fluid dynamics properties of the flow. CFD simulation is further utilized to visualize the flow regime in the radial direction.

The study showed that the Gnielinski and experimental Nusselt number trend in the same axial direction for inlet bulk temperatures further from the pseudo-critical temperature, however, cannot be used to effectively predict the experimental Nusselt number close to the pseudo-critical temperature. As the tube is inclined, the buoyant forces become negligible and there is no radial variation in temperature as the inclination angle approaches 90 deg. For downward flow, buoyant forces enhance the turbulent shear stress and heat transfer coefficient increases. It can be concluded that for inlet bulk temperatures near the pseudo-critical temperature, buoyant forces are stronger, but reduce as the inlet bulk temperature is increased. Maximum buoyant forces and heat transfer occur in the horizontal flow, and for inlet bulk temperatures near the pseudo-critical temperature. The results of this study provide insight that can be applied to design and development of pressure vessels, turbomachinery, and heat exchanger components in power plants or cycles that may use sCO2 as its working fluid.

Data Availability Statement

The datasets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request.

Nomenclature

Acs =

cross-sectional area of tube

Alocal =

local convection heat transfer area

Dext =

external diameter of tube = 12.57 mm

Dint =

internal diameter of tube = 10.92 mm

f =

Darcy friction factor

g =

acceleration due to gravity

Gr =

Grashof number

hb,i =

bulk enthalpy at station i

hb =

bulk enthalpy

HTClocal =

local heat transfer coefficient

I =

current flowing through test section

k =

thermal conductivity

kb =

thermal conductivity calculated using bulk CO2 properties

kSS =

thermal conductivity of Stainless steel

Llocal =

local length of heat transfer calculation domain

Ltotal =

total heated length of tube

m˙ =

mass flow rate

NuB =

local experimental Nusselt number at bottom surface

NuGn =

Gnielinski Nusselt number

Nulocal =

local experimental Nusselt number

NuT =

local experimental Nusselt number at top surface

Prb =

Prandtl number calculated using bulk CO2 properties

qvl =

volumetric heat generation in tube

q˙ =

heat flux

Qgen,electric =

electric power generated in test section

Qloss,radial =

heat loss to surrounding in radial direction

Reb =

Reynolds number calculated using bulk CO2 properties

Ri =

Richardson number

sCO2 =

supercritical Carbon Dioxide

Tamb =

ambient temperature

Tb =

bulk temperature

Tb,i =

bulk temperature at station i

Text =

temperature at external wall surface

Tint =

temperature at internal wall surface

Tw =

wall temperature

ui =

component of velocity vector in i direction

Uloss =

coefficient of radial heat loss

V =

voltage provided to test section

Vxy =

velocity magnitude in xy plane

xi =

cartesian coordinate in i direction

Δhb,i =

bulk enthalpy change at station i

ΔQCO2,local =

heat transferred to CO2 in heat transfer calculation domain

ΔTw −amb =

temperature difference between tube external wall and ambient

δi j =

Kronecker delta function

μ =

dynamic viscosity

μt =

turbulent viscosity

ρb =

bulk flow CO2 density

ρwall =

CO2 density at wall

ρ =

density

References

1.
Binotti
,
M.
,
Astolfi
,
M.
,
Campanari
,
S.
,
Manzolini
,
G.
, and
Silva
,
P.
,
2017
, “
Preliminary Assessment of sCO2 Cycles for Power Generation in CSP Solar Tower Plants
,”
Appl. Energy
,
204
, pp.
1007
1017
.10.1016/j.apenergy.2017.05.121
2.
Wu
,
C.
,
Xu
,
X.
,
Li
,
Q.
,
Li
,
J.
,
Wang
,
S.
, and
Liu
,
C.
,
2020
, “
Proposal and Assessment of a Combined Cooling and Power System Based on the Regenerative Supercritical Carbon Dioxide Brayton Cycle Integrated With an Absorption Refrigeration Cycle for Engine Waste Heat Recovery
,”
Energy Convers. Manage.
,
207
, p.
112527
.10.1016/j.enconman.2020.112527
3.
Ma
,
Y.
,
Liu
,
Z.
, and
Tian
,
H.
,
2013
, “
A Review of Transcritical Carbon Dioxide Heat Pump and Refrigeration Cycles
,”
Energy
,
55
, pp.
156
172
.10.1016/j.energy.2013.03.030
4.
Khadse
,
A.
,
Blanchette
,
L.
,
Kapat
,
J.
,
Vasu
,
S.
,
Hossain
,
J.
, and
Donazzolo
,
A.
,
2018
, “
Optimization of Supercritical CO2 Brayton Cycle for Simple Cycle Gas Turbines Exhaust Heat Recovery Using Genetic Algorithm
,”
ASME J. Energy Res. Technol.
,
140
(
7
), p.
071601
.10.1115/1.4039446
5.
Wright
,
S.
,
Davidson
,
C.
, and
Scammell
,
W.
,
2016
, “
Thermo-Economic Analysis of Four sCO2 Waste Heat Recovery Power Systems
,”
Fifth International sCO2 Symposium
,
San Antonio, TX
, Mar. 29, pp.
28
31
.
6.
Gabriel-Ohanu
,
E.
,
Khadse
,
A.
,
Vesely
,
L.
,
Raju
,
N.
,
Otto
,
M.
,
Kapat
,
J. S.
, and
Harris
,
K.
,
2021
, “
Optimization of a Primary Heat Exchanger for Flibe Molten Salt Nuclear Reactor With sCO2 Power System
,”
ASME
Paper No. GT2021-59939.10.1115/GT2021-59939
7.
Huai
,
X.
,
Koyama
,
S.
, and
Zhao
,
T.
,
2005
, “
An Experimental Study of Flow and Heat Transfer of Supercritical Carbon Dioxide in Multi-Port Mini Channels Under Cooling Conditions
,”
Chem. Eng. Sci.
,
60
(
12
), pp.
3337
3345
.10.1016/j.ces.2005.02.039
8.
Yang
,
C.
,
Xu
,
J.
,
Wang
,
X.
, and
Zhang
,
W.
,
2013
, “
Mixed Convective Flow and Heat Transfer of Supercritical CO2 in Circular Tubes at Various Inclination Angles
,”
Int. J. Heat Mass Transfer
,
64
, pp.
212
223
.10.1016/j.ijheatmasstransfer.2013.04.033
9.
Schnurr
,
N.
,
1969
, “
Heat Transfer to Carbon Dioxide in the Immediate Vicinity of the Critical Point
,”
ASME J. Heat Mass Transfer-Trans. ASME
,
91
(
1
), pp.
16
20
.10.1115/1.3580086
10.
Adebiyi
,
G. A.
, and
Hall
,
W.
,
1976
, “
Experimental Investigation of Heat Transfer to Supercritical Pressure Carbon Dioxide in a Horizontal Pipe
,”
Int. J. Heat Mass Transfer
,
19
(
7
), pp.
715
720
.10.1016/0017-9310(76)90123-X
11.
Kim
,
T. H.
,
Kwon
,
J. G.
,
Kim
,
M. H.
, and
Park
,
H. S.
,
2018
, “
Experimental Investigation on Validity of Buoyancy Parameters to Heat Transfer of CO2 at Supercritical Pressures in a Horizontal Tube
,”
Exp. Therm. Fluid Sci.
,
92
, pp.
222
230
.10.1016/j.expthermflusci.2017.11.024
12.
Bae
,
Y.-Y.
, and
Kim
,
H.-Y.
,
2009
, “
Convective Heat Transfer to CO2 at a Supercritical Pressure Flowing Vertically Upward in Tubes and an Annular Channel
,”
Exp. Therm. Fluid Sci.
,
33
(
2
), pp.
329
339
.10.1016/j.expthermflusci.2008.10.002
13.
Pidaparti
,
S. R.
,
McFarland
,
J. A.
,
Mikhaeil
,
M. M.
,
Anderson
,
M. H.
, and
Ranjan
,
D.
,
2015
, “
Investigation of Buoyancy Effects on Heat Transfer Characteristics of Supercritical Carbon Dioxide in Heating Mode
,”
ASME J. Nucl. Eng. Radiat. Sci.
,
1
(
3
), p.
031001
.10.1115/1.4029592
14.
Tanimizu
,
K.
, and
Sadr
,
R.
,
2016
, “
Experimental Investigation of Buoyancy Effects on Convection Heat Transfer of Supercritical CO2 Flow in a Horizontal Tube
,”
Heat Mass Transfer
,
52
(
4
), pp.
713
726
.10.1007/s00231-015-1580-9
15.
Guo
,
P.
,
Liu
,
S.
,
Yan
,
J.
,
Wang
,
J.
, and
Zhang
,
Q.
,
2020
, “
Experimental Study on Heat Transfer of Supercritical CO2 Flowing in a Mini Tube Under Heating Conditions
,”
Int. J. Heat Mass Transfer
,
153
, p.
119623
.10.1016/j.ijheatmasstransfer.2020.119623
16.
Liao
,
S.
, and
Zhao
,
T.
,
2002
, “
Measurements of Heat Transfer Coefficients From Supercritical Carbon Dioxide Flowing in Horizontal Mini/Micro Channels
,”
ASME J. Heat Mass Transfer-Trans. ASME
,
124
(
3
), pp.
413
420
.10.1115/1.1423906
17.
Jackson
,
J. D.
, and
Hall
,
W. B.
,
1979
, “
Forced Convection Heat Transfer to Fluids at Supercritical Pressure
,” Turbulent Forced Convection in Channels and Bundles: Theory and Applications to Heat Exchangers and Nuclear Reactors,
Hemisphere Publishing Corporation
,
New York
, pp.
563
612
.
18.
Lemmon
,
E.
,
Huber
,
M.
, and
McLinden
,
M.
,
2002
, “
Nist Standard Reference Database 23: Refprop Version 9.0
,”
National Institute of Standards and Technology
,
Boulder, CO
.
19.
Khadse
,
A.
,
Vesely
,
L.
,
Sherwood
,
J.
,
Curbelo
,
A.
,
Goyal
,
V.
,
Raju
,
N.
,
Kapat
,
J. S.
, and
Kim
,
W.
,
2020
, “
Study of Buoyancy Effects on Supercritical CO2 Heat Transfer in Circular Pipes
,”
ASME
Paper No. GT2020-15523.10.1115/GT2020-15523
20.
Franssen
,
J.-M.
, and
Real
,
P. V.
,
2012
, Fire Design of SteelStructures: Eurocode 1: Actions on Structures; Part 1-2: General Actions–Actions on Structures Exposed to Fire; Eurocode 3: Design of Steel Structures; Part 1-2: General Rules–Structural Fire Design,
John Wiley Ernst & Sohn
,
Berlin, Germany
.
21.
Bergman
,
T. L.
,
Incropera
,
F. P.
,
Lavine
,
A. S.
, and
DeWitt
,
D. P.
,
2011
, Introduction to Heat Transfer,
John Wiley & Sons
,
New York
.
22.
Kline
,
S. J.
,
1953
, “
Describing Uncertainty in Single Sample Experiments
,”
Mech. Eng.
,
75
, pp.
3
8
.
23.
Test Uncertainty
,
A.
,
2006
, “
PTC 19.1-2005
,”
Am. Soc. Mech. Eng.
,
3
, pp.
10016
5990
.
24.
Menter
,
F. R.
,
1994
, “
Two-Equation Eddy-Viscosity Turbulence Models for Engineering Applications
,”
AIAA J.
,
32
(
8
), pp.
1598
1605
.https://www.asme.org/getmedia/2534a386-a1b0-4059-89e7-7422a92c982f/35537.pdf
25.
Huang
,
D.
,
Wu
,
Z.
,
Sunden
,
B.
, and
Li
,
W.
,
2016
, “
A Brief Review on Convection Heat Transfer of Fluids at Supercritical Pressures in Tubes and the Recent Progress
,”
Appl. Energy
,
162
, pp.
494
505
.10.1016/j.apenergy.2015.10.080
26.
Lei
,
X.
,
Zhang
,
Q.
,
Zhang
,
J.
, and
Li
,
H.
,
2017
, “
Experimental and Numerical Investigation of Convective Heat Transfer of Supercritical Carbon Dioxide at Low Mass Fluxes
,”
Appl. Sci.
,
7
(
12
), p.
1260
.10.3390/app7121260
27.
Wang
,
X.
,
Xiang
,
M.
,
Huo
,
H.
, and
Liu
,
Q.
,
2018
, “
Numerical Study on Nonuniform Heat Transfer of Supercritical Pressure Carbon Dioxide During Cooling in Horizontal Circular Tube
,”
Appl. Therm. Eng.
,
141
, pp.
775
787
.10.1016/j.applthermaleng.2018.06.019