Numerical Simulation of Supercritical Carbon Dioxide Critical Flow in the Nozzle Tube

Supercritical fluids are characterized by a temperature and pressure higher than the critical state properties of these respective fluids. The supercritical carbon dioxide (S-CO₂) Brayton gas turbine cycle has been considered as an alternative fluid and has been used in numerous power cycle designs to increase the power density by reduction of system component sizes and increased thermal efficiency due to the higher operating temperature [1–3].

The coolant critical flow is treated as an important issue for nuclear reactor systems. In recent decades, numerous critical flow studies have been conducted for water under both two-phase and subcooled condition [2,4–8]. Supercritical carbon dioxide critical flow data about blowdown or rapid depressurization are favorable for certain aspects of safety analyses. However, experiments and theories for S-CO₂ are still insufficient. Mignot et al. [9] accomplished an experiment to measure critical mass flux of a long L/D tube in numerous stagnation thermodynamic conditions, geometry, and outlet tube roughness. The one-dimensional homogeneous equilibrium model was developed to predict critical flow. Edlebeck et al. [10] conducted an S-CO₂ flow test loop and measured supercritical and two phase CO₂ flow rate through short orifices. Orifices with diameters of 1 mm and orifice length-to-diameter ratios of 3.2 and 5 were tested. Now, some experiments are carried out to characterize S-CO₂ critical flow behavior; however, the analysis of the critical flow with computational fluid dynamics code is still required.

The critical flow in the nozzle tube is analyzed numerically with fluent 15.0. A mesh quality study is conducted for choosing the optimum meshing strategy. Experimental data [8–10] are used as a benchmark to examine the capability of the current approach. In addition, the effects of L/D and stagnation temperature (the temperature at which the fluid flow rate is zero) on critical mass flow are investigated. Also, the flow patterns inside the tube during critical flow are discussed based on simulation results.

Experimental data of Mignot et al. are used as a benchmark to examine the capability of the current approach. The test section can be divided into three parts, i.e., inlet part, tube part, and outlet part. Four test sections are analyzed. A two-dimensional, axisymmetric model is adopted in the calculation, as shown in Fig. 1. The geometrical size of the test sections is shown in Table 1. Pressure and temperature are given at the inlet boundary, and the pressure is given at the outlet boundary. The pressures and temperatures at the boundaries are set as the constant values shown in Table 2. The pressure at the outlet is set as 0.1 MPa. The transient flow field is calculated until a steady-state flow is observed.

For the computational fluid dynamics analysis, fluent ver.15.0 is used. The fluid is set as CO₂. The basic equations are the continuity equation, momentum equation, and energy equation. Because CO₂ in the nozzle is at high Reynolds number, the standard k-epsilon model and standard wall functions are selected to calculate turbulence. The pressure implicit with splitting of operators scheme is used as the flow solution method. The spatial gradient terms are discretized by the least-squares cell-based method. The pressure interpolation scheme is a second-order scheme. The convection terms of momentum and energy equations are discretized by the quadratic upwind interpolation for convective kinematics method for the higher accuracy calculation, and the convection terms of the other equations are discretized by the first-order upwind scheme. The simulation was iterated about 100,000 times to converge.

The Error₁ is the deviation between R–K real gas density and actual gas density. The Error₂ is the deviation between ideal gas density and actual density. By Eqs. (1)–(3), density is calculated, comparing their actual density with the density of the R–K real gas model. Table 3 shows the errors are less than 3.5%, confirming the R–K real gas model.

At first, results of a convergence study are presented to ensure the adequacy of the grid resolution. The test section no. 1 is selected for analysis. Five grid resolutions, 330, 925, 3700, 14,800, 59,200, are considered. A uniform structured grid strategy is used.

The stagnation pressure and temperature are set as 19.2 MPa, 263.7 °C. The outlet pressure is 0.1 MPa. The simulation results of different grid conditions are shown and compared with experiment result in Table 4.

Table 4 shows that the error decreases as the number of cells increases. As an optimum in computational time and error, a grid with 14,800 cells is used for the present calculation. The grid size is 0.32 mm × 0.32 mm.

The length-to-diameter ratio and inlet stagnation parameters are the main parameters affecting the critical flow. The simulation results of critical flow under different draw ratio and stagnation parameter conditions are shown in Table 5.

Table 5 shows that the error between simulation and experiment is less than 20%. Compared with experiments, the simulation result overestimates the critical mass flow rate. Figure 2 shows simulation and experimental results for the critical flow rate. As the length-to-diameter ratio increases, the mass flow rate decreases. For short test sections (L/D < 20), the critical mass flow rate decreases quickly with increasing length-to-diameter ratio. For long test sections (L/D > 20), the critical mass flow rate decreases slightly with increasing length-to-diameter ratio. The trend of simulation results of the length-to-diameter ratio effects agrees with the experimental data.

Some errors between experimental data and the simulation still exist, because of drastic parameter changes at inlet and outlet of the nozzle. The difference between the R–K real gas model and actual parameters also brings some errors to the results.

Table 6 shows the effect of stagnation temperature on critical flow when the L/D is 48.

Table 6 shows that the error between simulation and experiment is less than 25%. Compared with experiments, the simulation result overestimates the critical mass flow rate again. Figure 3 shows effect of stagnation temperature on critical flow rate when the stagnation pressure is 10.1 MPa and the L/D is 48. As stagnation temperature increases, the mass flow rate decreases. The predicted trend of the stagnation temperature effects agree with the experimental data. When the stagnation temperature increases, the fluid density decreases, which lead to less mass flow rate.

When the fluid flow enters the nozzle, there will be a flow area reduction at the inlet. The smallest cross section of the flow area is known as the “throat.” It results in an increase in the flow velocity and the velocity reaches a maximum in the throat.

In order to reveal the critical flow mechanism, the pressure and temperature distribution of case 9 is shown in Fig. 4 in detail.

At inlet, when the fluid is flowing through the tube, the pressure and temperature decrease quickly because of contraction and expansion. After tube inlet contraction, pressure and temperature increase a little due to flow area increase and velocity decreases lightly due to friction pressure losses.

For heat transfer, when critical flow happens, the velocity of CO₂ reaches local sonic speed in the nozzle. There is little time for heat transfer because of high speed. So, the temperature of the flow in the tube decreases slowly, as is shown in Fig. 4.

Case 4, case 9, and case 12 are selected for flow pattern behavior in this paper. The values of pressure and temperature along the central line of test section are presented in Fig. 5.

Three flow patterns in the nozzle tube during transient critical flow are obtained. From inlet to outlet of the tube, CO₂ may undergo follow phases in turn: (1) supercritical phase; (2) supercritical phase—gas phase; (3) supercritical phase—gas phase—liquid phase.

Results of case 9 and case 12 show that as the stagnation temperature decreases, the P-T path line is closer to the critical point and to the pseudocritical line. The simulation results of three cases show that, at outlet of the test section, CO₂ can be in liquid state and solid state due to low pressure and temperature. The temperature contour of case 9 is shown in Fig. 6. The temperature in the outlet of the test section changes dramatically. The flow behavior at the outlet of the test section in the experiment is shown in Fig. 4. The simulation result shows similar phenomena as the experiment (Fig. 7).

The critical flow in the nozzle tube is analyzed numerically based on fluent 15.0. Compared with experimental data, the simulation results overestimate the critical mass flow rate, and the error range is between 15% and 25%.

The simulation results show that, as L/D ratio increases, critical mass flow decreases. As the stagnation temperature increases, the critical mass flow decreases. Three flow patterns in the nozzle tube during transient critical flow are obtained and discussed. From inlet to outlet of the tube, CO₂ may undergo the following phases in turn: (1) supercritical phase; (2) supercritical phase—gas phase; (3) supercritical phase—gas phase—liquid phase.

This work is supported by National Natural Science Foundation of China (Grant Nos. 51506134 and 11605193).

a =

attractive potential among molecules, Pa m⁶ mol⁻²

b =

the R–K equation correction parameter introduced by the intrinsic volume of the molecule, m³/mol

L/D =

length-to-diameter ratio

P =

upstream pressure, Pa

P_C =

critical pressure, Pa

R =

molar gas constant, 8.31451 J/mol K

T =

upstream temperature, K

T_c =

critical temperature, K

V_m =

molar volume, m³/mol