Numerical Simulations of Cryogenic Hydrogen Cooling in Vortex Tubes with Smooth Transitions

Improving efficiency of hydrogen cooling in cryogenic conditions is important for the wider applications of hydrogen energy systems. The approach investigated in this study is based on a Ranque-Hilsch vortex tube (RHVT) that generates temperature separation in a working fluid. The simplicity of RHVT is also a valuable characteristic for cryogenic systems. In the present work, novel shapes of RHVT are computationally investigated with the goal to raise efficiency of the cooling process. Specifically, a smooth transition is arranged between a vortex chamber, where compressed gas is injected, and the main tube with two exit ports at the tube ends. Flow simulations have been carried out using STAR-CCM+ software with the real-gas Redlich-Kwong model for hydrogen at temperatures near 70 K. It is determined that a vortex tube with a smooth transition of moderate size manifests about 7% improvement of the cooling efficiency when compared vortex tubes that use traditional vortex chambers with stepped transitions and a no-chamber setup with direct gas injection.


Introduction
Hydrogen is often considered as the most promising candidate for the "green" fuel of the future. However, there are still hurdles that need to be overcome to make hydrogen sufficiently easy and economical to handle. In the liquid form, hydrogen becomes one of the most energy dense substances. However, in order to liquefy hydrogen, very low temperatures are needed, and this process is complicated by the fact that hydrogen undergoes exothermic transformation between ortho and para isomeric states during the cooling process. Therefore, improvement of the current and development of novel liquefaction systems for hydrogen is of high priority to the green energy community.
One of the promising mechanical components that can boost efficiency of cryogenic refrigeration systems is the so-called vortex tube. This device is known to produce two streams with significant temperature separation at two opposite outlets while taking a compressed fluid at the inlet (Figure 1a). Ranque [1] discovered this effect by injecting gas tangentially into a tube and producing a cold stream near the centerline while the hot gas was concentrated in the peripheral flow. Hilsch [2] investigated this phenomenon further, which led to the name Raque-Hilsch vortex tube commonly used nowadays for this device.
Being a mechanically simple system (no moving solid parts), vortex tubes have been extensively studied since the 1950s, and several reviews can be found in [3][4][5]. A number of studies addressed physical mechanisms for temperature separation, characterization of the system performance and some optimization. Among various explanations for the vortex-tube thermal effect, the most common include shear work done by shear work of the inner vortex on the outer flow exposed to wall friction [6], gas expansion near the tube centerline due to rotation [7] and turbulent fluctuations in the radial direction in However, there is still room for improvement of the vortex tube configurations. For hydrogen liquefaction processes, even small increases in efficiency may contribute to significant benefit, when such systems are implemented on a mass scale. An initial effort toward the development of vortex tubes for cryogenic hydrogen systems was discussed in [11].
The present study focuses on one novel shape modification of the vortex tube associated with the smooth transition between the vortex chamber and the main straight tube (Figure 1b). It was discovered that rounding of corners at this transition [12,13], making the entire tube convergent [14,15] and using conical shapes [16] can lead to the improvement of air-based vortex tubes. Experimental optimization of structural parameters of a control valve and a hot tube in a vortex-tube air separator and numerical studies of internal flow features were conducted in [17]. Performance of a novel convergent vortex tube was parametrically investigated in a test series [18], where it was found that much shorter tube lengths are possible in convergent setups in comparison with straight tubes and that valves can be eliminated as well. In the present work, the effect of a smooth transition between the vortex chamber and the main tube in a system using gaseous hydrogen as the working fluid is investigated with help of computational fluid dynamics (CFD) simulations.
In one of the initial CFD studies with vortex tubes [19], a simplistic axisymmetric flow with swirl was modeled, and some adjustment of geometrical parameters in CFD was needed to achieve agreement with test data. Computational studies using economical Reynolds-averaged Navier-Stokes (RANS) models generally suggested that the standard k − ε turbulence model gives better prediction [16,20]. More expensive LES-based simulations provided more insight on flow patterns [21,22].
Applications of CFD programs for modeling fluid phenomena in hydrogen and cryogenic systems have steadily increased in the last two decades. Several papers on numerical simulations relevant to cryogenic vortex tubes have been published. For example, energy and species separation in air at cryogenic conditions were investigated in [23]. Hydrogen-filled vortex tubes with traditional geometry were modeled in [24] for a wide range of inlet temperatures, including cryogenic states; and it was found that the vortex tube performance generally decreases at low temperatures. In the current work, highfidelity numerical simulations are employed to simulate cooling of cryogenic gaseous hydrogen in a vortex tube with novel geometrical modifications. A short version of the present paper has been previously reported at the conference [25].

Computational Aspects
Modeling of vortex tubes in this study was carried out with the CFD program STAR-CCM+ (Siemens, Munich, Germany), version 13.04.011. It employs a finite-volume coupled flow solver with the second order discretization. For computational economy, only steady solutions were sought, and simplified turbulence models based on the RANS formulation were utilized. The governing Reynolds-averaged equations include the mass, momentum and energy equations, where u i are the time-averaged velocity components, p is the pressure, ρ is the fluid density, µ is the fluid viscosity, −ρu ı u  is the turbulent stress tensor, E is the total energy per unit mass and q j are the heat fluxes.
To model the turbulent stress, several RANS turbulent models were investigated, and the standard k − ε model was found to perform best in the present study, as discussed in the next section. Moreover, the two-layer, all Y+ option of this model was employed due to significant flow variability in the vortex tube. The specific implementation details of the numerical methods in the applied software can be found in STAR-CCM+ Manual [26].
Since the current study analyzes flow of gaseous hydrogen at cryogenic conditions, a real-gas equation of state was employed to describe the gas properties. Among several approaches available in the literature for real-gas modeling, the Redlich-Kwong model was found to treat equilibrium hydrogen gas quite accurately in the range of properties of interest to this study [27,28]. Hence, the Redlich-Kwong model was adopted for all simulations shown in this paper.

Model Verification and Validation
The verification study has been conducted with the most promising configuration among considered in this paper. The tube geometry with inlets and outlets is illustrated in Figure 2, while its dimensions are given in the next section. There are four inlets where stagnation properties are specified (3 atm and 77 K), the cold outlet has the stagnation pressure of 1 atm, and the hot-outlet flow rate is selected to make the cold-flow fraction equal to ζ = 0.46. This fraction is defined as the ratio of the mass flow rates at the cold outlet, . m C , and at the inlet, The other domain boundaries are solid walls which are treated as adiabatic no-slip surfaces.
The unstructured numerical grid was built inside the domain containing polyhedral cells that work well in situations with complex flows. Near the solid walls, five prismatic cell layers were constructed to cover the boundary layer. Because of high variations of the flow velocity inside vortex tubes and associated difficulty of ensuring the same cell thickness in terms of Y+ values, the two-layer, all-Y+ option was activated in all turbulent models. This approach allows for resolving a viscous sub-layer when the cell thickness corresponds to Y+ below 5, while relying on wall functions for large Y+ values and using a combined model for intermediate Y+ values.
Three numerical grids (fine, intermediate, coarse) were constructed by scaling the cell linear dimension with a factor of 2 between neighboring grids. An illustration of the fine-grid slices in several tube sections is given in Figure 3. The cell count on the fine grid was about 1.6 million.  In the grid independence study, the difference between total temperatures at the inlet and cold outlet was used as the metric of convergence. The obtained results for T in,tot − T C,tot at three grids are given in Figure 4. One can see that this parameter approaches a constant value as the number of cells increases. The numerical uncertainty, evaluated with help of Richardson extrapolation [29] and factors of safety [30], produced an estimated error of about 10% relative to results obtained with the fine grid. In all following simulations carried out in this study, the same fine-mesh settings were employed when constructing numerical grids. In this work, a validation study was run against test data available for vortex tubes with air at normal conditions [31]. They reported the temperature differences between the hot and cold outlets for a common vortex tube configuration of relatively low performance. The comparison between the present numerical results and previous experimental data are shown in Figure 5. As once can see, simulation results obtained with the standard k − ε turbulence model are in satisfactory agreement with experimental data, while the realizable k − ε model predicts much lower temperature differentials. Other turbulence models were also tried, including the standard and (shear stress transport) SST k − ω models and Reynolds stress transport (RST) model. They produced results in between the two numerical models shown in Figure 5. As the standard k − ε turbulence model performed best with other numerical settings in the present study, it was adopted for all other simulations.

Simulation Results
To come up with a reasonably efficient initial design of a vortex tube, the semiempirical procedure of Merkulov [3] was followed. This method involves experimental charts and correlations that were obtained from tests of various air-filled vortex tubes. Here, the gas properties for hydrogen at 77 K were used instead. The temperature of 77 K corresponds to the boiling point of liquid nitrogen, which allows for relatively inexpensive cooling of hydrogen down to this temperature. Further cooling would require expensive refrigerants, such as helium-neon mixtures, and that is where the vortex tube application becomes attractive.
The external input parameters were selected as follows: (i) the inlet stagnation pressure and temperature are 3 atm and 77 K, (ii) the cold outlet pressure is 1 atm, (iii) the desirable cold-outlet stagnation temperature is about 70 K (since the expected temperature drop in vortex tubes with pressure ratio of 3 is about 10%) and (iv) the desired input mass flow rate is around 5 g/s. The applied design method produced the following recommendations for the device geometry: the main tube length of 133 mm and diameter of 14.9 mm, the cold-outlet diameter of 7.3 mm, and the total inlet area of 16 mm 2 . To ensure adequate performance at this relatively low length-to-diameter ratio, vortex stopper (de-twister) is needed near the hot outlet, which was shown in Figure 2. Since the vortex-tube performance usually benefits from several distributed inlets, four square nozzles of side 2 mm have been chosen here. The optimal cold-flow fraction suggested by this method for the selected geometry is about 0.46.
The initial tube configuration was selected to have a constant-diameter vortex tube with direct tangential nozzle integration at one end of the main tube. This setup is illustrated as case A in Figure 6. Since a stepped vortex chamber is commonly used in practice, such a system was created by displacing nozzles outward from the tube centerline by two nozzle side lengths (case B in Figure 6). The next configuration with a smooth transition between the vortex chamber and the main tube was considered, as shown by case C in Figure 6. (It was also used for the verification study discussed above.) Since this geometry demonstrated performance improvement, even wider vortex chamber was also tried with diameter greater than the main tube diameter by 8 mm (four nozzle sides), as shown by case D in Figure 6. These four configurations were meshed with a fine mesh following the same approach as in the above validation study. Simulations were then carried out for these cases at the same total pressure at the inlet and cold outlet and the same inlet total temperature. The hot-outlet pressure was adjusted to keep the cold-flow fraction at the same value (0.46).
The main dimensional characteristics of the vortex tubes include the stagnation temperature separations at the outlets in comparison with the inlet stagnation temperature, T in,tot − T C,tot and T H,tot − T in,tot . Another common metric for the vortex-tube performance is the isentropic temperature efficiency that compares actual temperature drop at the cold end to that in the isentropic processes at the same pressure ratio p r = P in,tot /P C,tot .
These characteristics obtained from simulations with four tube configurations are given in Figure 7. It is apparent that introducing a stepped vortex chamber (case B) reduces the system efficiency in comparison with a straight tube (case B). However, employing a smooth transition between the shifted nozzles and the main tube (case C) helps boost the temperature separation (up to 9 K at the cold side) and efficiency (η T ≈ 0.44). Case C seems to have expansion of the tube near the inlets closest to the optimum (among the studied cases), since case D with even larger cross-sectional area exhibits performance degradation. Some flow characteristics of the four studied cases are presented in Figures 8-10. The flow patterns shown with help of line convolution integrals, as well as velocity magnitudes, are given in Figure 8 at longitudinal tube sections. The maximum velocity in the predominantly tangential direction is observed near the inlets. The incoming flow forms a swirl that propagates downstream keeping highest velocities closer to the walls (although velocity at the wall surface drops to zero due to no-slip condition). The flow regions near the centerline are relatively slow. The lowest velocities are shown near the hot outlet where the vortex stopper significantly decelerates the flow. It can be noticed that the vortex chamber expansions alter the size of the high-velocity swirl. It appears that case C with smooth transition and moderate widening keeps the fast near-wall flow the longest distance along the tube, which eventually leads to greater temperature separation.  In Figure 8, one can also notice slices of quasi-toroidal vortices between the outer flow with axial direction toward the hot end and the inner flow with axial motion to the cold outlet. In cases A and B these vortices are rather thin and elongated. With wider tube enlargements, these vortices grow in the radial direction, which may help to keep the faster outer flow in case C. In Case D, several vortex systems are present, which likely leads to excessive dissipation and reduction in the system efficiency. The flow Mach numbers in the tube transverse sections passing through the middle of nozzles are depicted in Figure 9. The highest Mach numbers (around one) are observed in the original tube (case A) without a vortex chamber. These near-sonic zones are formed at short distances from the entrances where the incoming flow from the nearest nozzle is merging with the flow supplied though other nozzles. Near the tube centerlines, the flow velocity drops substantially. In case B with a stepped vortex chamber, the high-Mach number regions are shifted toward the tube center, where the flow from the inlets enters the main tube. More uniform distribution of Mach numbers is observed in case C, while their magnitudes still reach high values. Even more uniform, but also substantially lower Mach numbers, appear in the widest case D.
The selected streamlines with superimposed total temperature magnitudes are illustrated for four studied cases in Figure 10. All streamlines have swirl shapes due to high tangential velocities. The spatial period in the axial directions increase away from the nozzles indicating reduction of the tangential speeds, while the axial velocity changes little. The streamlines in the outer region generally demonstrate an increase in T tot as they approach the hot outlet, with the exception of case D where the total temperature peaks at the end of the smooth transition between the inlets and the main tube. The streamlines originated at the nozzle closer to the tube centerline initially propagate toward the hot end, but then they change direction and proceed toward the cold outlet while reducing their radii.
Additional simulations have been conducted for cases A (original) and C (highest performing among studied setups) at low and high cold-flow fractions. The results for the cold temperature separation and isentropic efficiency are shown in Figure 11. Besides these parameters, another metric important in practical applications is introduced by multiplying the cold-flow fraction ζ and the isentropic efficiency η T . This product characterizes the cooling capacity of the vortex tube rather than just the cold temperature. As one can see in Figure 11, case C slightly outperforms the original setup A in the entire range studied. The isentropic efficiency is the highest at the low ζ, while the capacity metric peaks at the high cold-flow fraction. Figure 11. Metrics of performance of vortex tubes: (a) temperature separation at the cold outlet, (b) isentropic temperature efficiency, (c) non-dimensional cooling capacity. Stars, highest performing (among studied setups) case C; triangles, original case A.

Concluding Remarks
In the present study, flow of cryogenic hydrogen gas inside vortex tubes was numerically simulated. Verification and validation indicate reasonably good accuracy of this modeling. The innovative geometry of the vortex tube, involving a smooth transition between the vortex chamber and the main tube, resulted in increased by about 7% performance of such a tube relative to a standard straight tube. This improvement can be associated with slower degradation of the vortex strength along the tube away from the inlets, which enhances the energy transfer from the core to peripheral fluid regions. However, if the tube widening near the entrance becomes excessive, the efficiency drops since the inlet flow becomes insufficient to maintain high-velocity swirl in a wide tube.
Further optimization of the vortex tube may involve variations of other geometrical parameters of the system. A next step in modeling vortex tubes intended specifically for hydrogen liquefaction can include even lower temperatures and possibly involve hydrogen condensation and transition between ortho/para phases of hydrogen. Hence, actual properties of hydrogen beyond the simplified equations of state will need to be used. Modeling of condensation may require finer numerical grids and multi-phase approaches.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.

Conflicts of Interest:
The authors declare no conflict of interest.