1. Introduction
As one of the main development directions in modern internal combustion engines, high-speed and high-power density internal combustion engines are characterized by high injection pressure and speed and short injection duration. With the increasing injection pressure in a common high-pressure rail system, the cavitation problem inside the injector nozzle holes becomes more and more serious. When the pressure difference ΔP between the inlet and outlet of the injector hole is large, a low-pressure area will be generated near the wall of the flow path inside the injector hole, resulting in the pressure of high-speed fuel in this area rapidly dropping to its saturated vapor pressure value, thus leading to cavitation. Cavitation in the injector hole affects the subsequent spray combustion effect, and long-term cavitation corrosion will lead to spray hole wall peeling, affecting the injection accuracy and reliability of engine operation. Currently, the use of computational fluid dynamics (CFD) methods to validate relevant injector fuel cavitation flow tests and thus study the cavitation in injector holes is one of the hot research topics in the field of internal combustion engines.
The current fuel cavitation tests that are widely used to verify the reliability of numerical model calculations are mainly the two-dimensional injector hole test proposed by Sou et al. [
1] and the microchannel test proposed by Winklhofer et al. [
2,
3]. Because the microchannel test is closer to the actual work of the injector hole, in recent years, scholars at home and abroad have widely used it to validate various new numerical models of the spray hole. Dai et al. [
4] used the Winklhofer test for the validation of the reliability of the model of a curved structure nozzle and used the nozzle model to study the effect of nozzle curvature on the transient cavitation characteristics of the nozzle, finding that the nozzle curvature can significantly reduce the cavitation in the hole degree. The greater the arc, the smaller the area of cavitation, and the smaller the average velocity of the spray hole exit, but the mass flow rate did not change significantly. By combining the Winklhofer test data, Sa et al. [
5] used the realizable
k-ε turbulence model and Schnerr and Sauer (SS) cavitation model in their calculations and verified that the cavitation distribution of the constructed multiphase flow model was consistent with the test, and then quantitatively investigated the effect of the spiral counter-groove structure on the injector needle valve on the turbulent fuel flow in the orifice and the subsequent fuel injection. Cristofaro et al. [
6] used the large eddy simulation (LES) model to calculate the cavitation distribution of his multiphase flow model, in agreement with the experimental data, and investigated the effect of liquid fuel viscosity on the mass flow rate, velocity profile and cavitation distribution during the microchannel throttling flow. Zhao et al. [
7] used a model with a grid size of 10 μm to verify the cavitation distribution, outlet mass flow rate and cross-sectional velocity distribution of the model similarly to the experiment, and finally investigated the effect of the compressibility of the fuel on the emission coefficient and critical cavitation number of the fuel in the nozzle. Guo et al. [
8] used the homogeneous relaxation model and the RNG
k-ε turbulence model calculations to verify that the model used was consistent with the experimental outlet mass flow rate and that the error in the cavitation distribution was minimal, thus investigating the effect of needle valve motion on the development of cavitation in the nozzle.
Although numerical calculation results closer to the experimental phenomena can be obtained in the study of the Winklhofer test at present, the influence of the selected computational model on the computational results and experimental validation has not been reported in the above-mentioned papers because the computational models and meshing strategies used are not the same. In recent years, some scholars have conducted related studies based on the influence of cavitation model selection on simulation results. Saha et al. [
9] investigated the influence of three cavitation models such as Zwart-Gerber-Belamri (ZGB), SS and Saha-Abu-Ramadan-Li (SAL) on simulation results with the selection of single-fluid and two-fluid methods. Piehl et al. [
10] carried out similar research work using three cavitation models such as Original Fire and SS, and Singhal et al. and Kumar et al. [
11] combined multiple URANS models with ZGB models to investigate the effect of the selection of Unsteady Reynolds Averaged Navier-Stokes (URANS) models on the computational results and verified them by comparing with the experimental data of Duke et al. [
12]. He et al. [
13] evaluated the predictive ability of different models for string cavitation in nozzles based on turbulence models such as RNG
k-ε, SST
k-ω, WALE LES and VLES, and found that LES and VLES obtained more accurate simulation results when studying the above problem. However, with the gradual use of LES and DES in microchannel cavitation simulations, it became more important to compare the differences in the effects of LES, DES, and URANS models on microchannel cavitation flow. The Winklhofer microchannel test mainly has four result indicators: pressure gradient distribution, cavitation distribution, outlet mass flow rate and cross-sectional flow rate distribution. In order to obtain more reliable research conclusions, the simulation data need to be compared and verified comprehensively with the above-proposed test parameters.
In this paper, we used Fluent 2020R2 commercial software to numerically simulate the pressure gradient, cavitation distribution, outlet mass flow rate and cross-sectional flow velocity distribution under various differential pressure conditions in the Winklhofer microchannel test based on five combinations of turbulence and cavitation models in turn, to compare and discuss the differences in numerical results among the model combinations and to evaluate and summarize the influence of model selection on the test validation results.
4. Calculation Results and Discussion
When fuel flows through the microchannel, the flow channel cross-section suddenly becomes smaller, the flow velocity increases and the pressure drops below the saturated vapor pressure of fuel, at which time flow separation and cavitation will occur. These phenomena directly affect the pressure and velocity distribution in the microchannel. They are reflected in the resultant parameters such as outlet mass flow rate, cavitation area distribution, flow velocity distribution and pressure gradient distribution. In this paper, the distribution characteristics of the above physical parameters are studied separately, the differences in results under multiple numerical model application scenarios are compared and analyzed, and finally, the accuracy of the numerical calculation of compressible fuel under multiple model applications is discussed.
4.1. Comparison of Pressure Gradient
In order to elucidate the intrinsic mechanism of cavitation incipient to cavitation development in the channel, this paper first analyzes the pressure gradient variation along the flow direction at the centerline of the model under the ΔP = 70 bar, as shown in
Figure 5. Under this pressure condition, the fuel flow in the microchannel reaches the critical cavitation condition. However, the fuel pressure value first decreases rapidly from the channel inlet to the middle of the cavitation region because the cavitation region is not fully developed at this time. It drops to a small value in the middle of the cavitation region with the smallest local basin cross-section, and then the basin cross-section increases and the pressure value rises back. After the fuel fluid flows through the end of the cavitation region and reattaches back to the wall, the basin cross-section gradually decreases and the flow rate gradually increases as the fluid continues to approach the channel exit, resulting in a gradual decrease in pressure as well.
Figure 5 illustrates the differences in pressure gradient profile changes for multiple turbulence and cavitation model applications. The pressure gradient curves obtained from different model combinations almost overlap and are similar to the experimental results as the fuel flows from the upstream buffer basin to the vicinity of the microchannel inlet. As the fuel flows through the cavitation and non-cavitation regions, all five model combinations, except for the RNG
k-ε and ZGB model combinations, show a trend where the pressure value first drops to a minimal value, then rises, and finally decreases. The errors between the pressure gradient results of the five models and the experimental data in each zone are given in
Table 3. As can be seen from
Table 3, except for the model combinations of LES and ZGB, the results obtained in the two zones of −1.00~−0.25 mm and −0.25~0.22 mm are close to the experimental data, and the errors are within 6%. In the 0.22–0.50 mm range, the best simulation results were obtained for the RNG
k-ε and ZGB model combination with an error within 6%; in the 0.50–2.00 mm range, the simulation results for the RNG
k-ε and ZGB, DES and ZGB model combinations were within 4% error. In general, the error between the pressure gradient simulation results obtained from the combination of RNG
k-ε and ZGB models and the experimental data is the smallest. From Equation (9), it can be seen that the RNG
k-ε model is based on the Standard
k-ε model with the addition of terms that can improve the accuracy for fast strain flow [
16]. Therefore, the RNG
k-ε model shows superior simulation performance in the two main zones of fast strain flow in microchannels, −0.25 to 0.22 mm and 0.22 to 0.50 mm.
4.2. Comparison of Cavitation Distribution
Figure 6 compares the cavitation distribution results between the experimental results and the simulation of various model combinations. Due to the ΔP between the front and rear buffer regions and the increase in fuel flow rate in the microchannel, the fuel pressure in the local region decreases to the fuel evaporation pressure, and the fuel in the region instantly evaporates into the gas phase, and the cavitation occurs.
From the vapor volume fraction cloud diagram, it can be seen that there are significant differences between the cloud diagrams obtained from different model combinations and the experimental results. Due to the presence of non-condensable gases in the actual fuel flow, the fuel liquid is more prone to gas nucleation during the pressure drop, and cavitation is more likely to occur. At 60 bar ΔP, the cavitation at the channel’s entrance was simulated by all five model combinations, consistent with the experimental results. For 70 bar ΔP, the cavitation clouds obtained by the LES and ZGB model combinations were close to the experimental observations. In contrast, the cavitation domain lengths obtained by the other five model combinations were smaller than those in the experimental results. For the pressure difference of 80 bar, the simulation results of the LES and ZGB and DES and ZGB model combinations were able to better reflect the experimental data, while the cavitation domain lengths obtained from the SST - and ZGB model combinations were slightly smaller than the experimental values.
It can be seen from the experimental results that as the ΔP increased from 60 to 80 bar, the gas phase fuel was no longer confined to the reflux region at the entrance of the microchannel, but gradually extended along the wall to the exit. However, the cavitation domains observed in the microchannels in the experiments were still larger than those obtained in the simulations, and the cavitation distributions obtained from different turbulence models were different. Yu et al. [
22] showed that different cavitation distributions were observed using different cavitation models for the same microchannel structure. Altimira et al. [
23] also found that the cavitation domain was larger in the experiment than in the simulation and pointed out that fuel compressibility affects the vapor condensation rate, which may be the reason for the small cavitation domain in the simulation. Meanwhile, he further found that using the different densities of grid distribution and adjusting the size of control parameters affecting cavitation bubble development in the cavitation model did not significantly affect the cavitation distribution at the near-wall surface. Giannadakis et al. [
24] showed that the cavitation domain development did not depend on the initial diameter of the cavitation bubble.
In summary, the cavitation domain attached to the microchannel wall developed as the ΔP between inlet and outlet increased, gradually extending along the wall to the outlet of the channel. The reason for the variability in the cavitation distribution between different model combinations at the same inlet and outlet pressure difference was the difference in the mathematical formulation of the different turbulence models because the cavitation domain was mainly composed of multiple large vortices clustered at the macroscopic level (one can see the wake of the vortices attached from the wall in the cavitation cloud at 60 bar pressure difference from the Winklhofer experiment). The above error analysis process also reflects that the LES model directly resolves the large vortices in a more suitable way for the simulation calculation of the cavitation domain than the way the remaining four turbulence models deal with the flow.
4.3. Comparison of Outlet Mass Flow Rates
Figure 7 reveals the microchannel outlet mass flow rate’s variation characteristics and the differences in the simulation results obtained from five model combinations. In this paper, the outlet mass flow rate is studied for thirteen ΔP conditions (ΔP
1 = 19 bar, ΔP
2 = 45 bar, ΔP
3 = 58 bar, ΔP
4 = 60 bar, ΔP
5 = 63 bar, ΔP
6 = 65 bar, ΔP
7 = 67 bar, ΔP
8 = 69 bar, ΔP
9 = 70 bar, ΔP
10 = 71 bar, ΔP
11 = 75 bar, ΔP
12 = 80 bar, ΔP
13 = 85 bar), which are entirely consistent with the conditions used in the experiments.
From the experimental results, it can be seen that in the process of gradually increasing the ΔP from 19 bar to 85 bar, the outlet mass flow rate firstly increases with the increase in ΔP between the inlet and outlet. Then, the mass flow rate curve shows an inflection point at 70 bar ΔP, after which the outlet mass flow rate maintains a certain value and no longer changes with the increase in ΔP. As shown in
Figure 7, numerical results obtained by all five model combinations were in good agreement with the experimental data in the range of 19 bar to 70 bar ΔP. As can also be seen from
Table 4, by dividing the model results into two zones (19 to 70 bar and 70 to 85 bar) with different variation trends of outlet mass flow rate according to the experimental data, the error between the results obtained from the five models and the experiments were within 6%. Among them, the simulation error of the LES and ZGB model combination was the smallest, and the error was within 2% in both zones. After a pressure difference of 70 bar, the results of the four model combinations, except for the SST
-
and ZGB model combinations, were slightly larger than the experimental data, but the trend of no further increase in the outlet mass flow rate could be reproduced. The LES (WALE) model assumes zero turbulent viscosity for the laminar shear flow, an assumption that allows the model to treat the laminar region in the domain correctly [
18]. As a result, the LES(WALE) model yielded simulation results that were in general agreement with the experimental data during the process in which the microchannel undergoes cavitation priming at 60 bar ΔP and subsequently progresses to sheet cavitation at 70 bar ΔP.
4.4. Comparison of Velocity Distribution
The velocity distribution along the
y-axis for two cross-sections (x
1 = 0.053 mm; x
2 = 0.17 mm) within the microchannel at 55 bar and 67 bar ΔP are shown in
Figure 8. Since the gasified fuel occupies the cavitation region when the ΔP is 67 bar, it is not easy to measure the velocity variation at the near-wall surface during the experiment [
5]. According to the trend of cross-sectional flow velocity distribution in the Winklhofer test, the simulation results of cross-sections x
1 and x
2 under two pressure differences can be divided into sections to calculate the average value of flow velocity and tabulated and summarized, as shown in
Table 5,
Table 6,
Table 7 and
Table 8.
Under the 55 bar ΔP condition, the simulation results at x
1 can somewhat reflect the velocity variation trend in the experimental data, as shown in
Figure 8b,d. On the other hand, at x
2, all model combinations cannot reproduce the complex flow velocity variation in the range of about 50~100 μm away from the upper and lower walls of the microchannel. As can be seen from
Table 5 and
Table 7, the simulation results are still quite different from the experimental data in the near-wall basin at a distance of about 0~50 μm from the upper and lower walls. When the
y-axis position point is about 50~250 μm (i.e., the central basin of cross-sections x
1 and x
2), the results obtained from the RNG
-
and ZGB models are closest to the experimental data, i.e., the calculation errors in both cross-sections are within 11%. The calculation errors of the DES and ZGB models are the largest, and the calculation errors in both cross-sections are over 14%.
Combined with
Figure 8c,e, the results obtained from the combination of the five models all match the trend of the test data at section x
1 under the 67 bar ΔP condition, but this still makes it difficult to reproduce the flow velocity distribution in the range of about 50–100 μm away from the wall at x
2. As can be seen from
Table 6 and
Table 8, the computational error of the flow velocity distribution near the wall at the 67 bar ΔP condition is slightly improved compared to that at the 55 bar ΔP condition. However, the simulation results are still somewhat different from the experimental data. For the central basin of cross-sections x
1 and x
2, the calculation errors of RNG
-
and ZGB models are still the smallest, i.e., the calculation errors at both cross-sections are within 8%. In contrast, DES and ZGB have the largest computational errors, which exceed 19% in both cross-sections.
Since the RNG
-
model can construct an effective viscosity term based on the RNG theory that takes into account the low Reynolds number effect [
16], this means that the RNG
-
model has superior performance in dealing with the near-wall flow or microchannel flow where the wall surface has a strong influence on the flow. Furthermore, the RNG
-
model retains the computational accuracy of the
-
type model for the high Reynolds number flow region. Therefore, we believe that the model treatment described above is the main reason why the RNG
-
and ZGB models eventually yield the simulation results closest to the experimental data.