# Modeling and Simulation of Diesel Injection at Transcritical Conditions

## 1 Introduction

The thermodynamic analysis of diesel engines shows that a higher chamber pressure will yield not only higher energy output but also higher engine efficiency. The need for improved engine efficiencies has motivated substantial research interests in high-pressure combustion systems, in which operating conditions achieve and exceed critical conditions. In diesel engines, typically a subcritical liquid fuel is injected into a supercritical ambient gaseous environment, and the fuel jet goes breakup, heating, and mixing processes before combustion. This process is referred to as a transcritical injection, which is schematically shown in \creffig:subsupercritical in comparison with a classical subcritical injection process. At elevated pressure, the fluid properties exhibit liquid-like densities and gas-like diffusivities, and the surface tension and enthalpy of vaporization approach zero [1]. This phenomenon has been shown by recent experimental works [2, 3, 4] and theoretical analysis [5, 6]. However, these complex processes are still not fully understood. Therefore, to provide insights into the high-pressure combustion systems, accurate and robust simulation tools are required for the characterization of supercritical and transcritical flows where large density gradients and thermodynamic anomalies occur.

The transcritical nature of flows in a diesel injection process has motivated several recent studies to utilize the diffused interface method for the modeling of the injection sequence. In contrast to a sharp interface method, such as a volume of fluid method, where interfaces are explicitly tracked or resolved in the computation domain, the diffused interface method artificially diffuses the interfaces. This is particularly attractive for transcritical flows where interfaces are not present. However, it remains an open research question whether interfacial flows or droplets exist under conditions relevant to real applications [1, 7]. Lacaze et al. [8] performed a large-eddy simulation (LES) using a real-fluid state equation to analyze the details of the transient mixing and the processes leading to autoignition during diesel injection. It was found that the large density ratio leads to significant penetration of the jet and enhanced entrainment effects. Ma et al. [9] developed a numerical scheme for the robust simulation under diesel relevant conditions. Knudsen et al. [10] employed a compressible Eulerian numerical model to describe the liquid fuel injection process under consideration of the internal nozzle flow and compressibility effects. Matheis and Hickel [11] developed a thermodynamic model where phase separation is considered through vapor liquid equilibrium calculations to enable the robustness of the diffused interface method when a fully conservative scheme is used.

Several challenges require considerations for the modeling of transcritical flows using the diffused interface method. One is the large density gradients caused by the liquid-like density in the dense fluid region and the gas-like density in the supercritical ambient gas region. A typical treatment is to introduce numerical dissipation locally where large density gradients occur [12, 13, 14, 15, 16]. Another difficulty arises from spurious pressure oscillations that are generated from nonlinearities of the real-fluid equation of state. These oscillations arise when fully conservative schemes are applied to multicomponent flows [17]. Schmitt et al. [13, 18] derived a quasi-conservative scheme by relating the artificial dissipation in the conservation equations and setting the pressure differential to zero. Terashima et al. [14, 19] solved the pressure equation instead of the energy equation. Ma et al. [16] developed a double-flux model with entropy-stability to eliminate spurious pressure oscillations in transcritical flows.

In this study, a recently developed compressible real-fluid solver [9, 16] is used for the simulation of a diesel injection process under realistic engine conditions. The governing equations of the conservation laws will be presented followed by the thermodynamic and transport models used for the description of the transcritical flows. The numerical methods employed in this study will be briefly discussed. These simulation results will be compared to available experimental measurements to examine the performance of the current numerical scheme and the diffused interface method for the prediction of the flow structures and mixing behaviors in a diesel injection process. The paper finishes with conclusions.

## 2 Mathematical Model

### 2.1 Governing Equations

The governing equations are the conservation laws for mass, momentum, total energy, and species, taking the following form:

(1a) | ||||

(1b) | ||||

(1c) | ||||

(1d) |

where is the density, is the velocity vector, is the pressure, is the specific total energy, is the mass fraction of species , is the diffusion coefficient for species , and in which is the number of species. The viscous stress tensor and heat flux are written as:

(2a) | ||||

(2b) |

where is the temperature, is the dynamic viscosity, is the thermal conductivity, and is the partial enthalpy of species . The total energy is related to the internal energy and the kinetic energy:

(3) |

The system is closed with a state equation, which is here written in the general form as,

(4) |

The formulation of the state equation are discussed next.

### 2.2 Equation of State

For computational efficiency and the accurate representation near the critical point [20], the Peng-Robinson cubic state equation [21, 22] is used in this study, which can be written as:

(5) |

where is the gas constant, is the specific volume, and the coefficients and are dependent on temperature and composition to account for effects of intermolecular forces.

The extended corresponding-states principle and pure fluid assumption for mixtures are adopted in this study [23, 24]. Mixing rules are applied for coefficients and in \crefEQ_PR_EQUATION_P:

(6a) | ||||

(6b) |

where is the mole fraction of species . The coefficients and are evaluated using the recommended mixing rules by Harstad et al. [25]:

(7a) | ||||

(7b) | ||||

(7c) |

The critical mixture temperature, pressure, molar volume, compressibility, and acentric factor are determined according to Harstad et al. [25].

### 2.3 Transport Properties

The dynamic viscosity and thermal conductivity are evaluated using Chung’s method with high-pressure correction [27, 28]. However, the evaluation of the viscosity in this model in known to produce oscillations for certain mixtures [15, 29]. To resolve this issue, in the present study the acentric factor is set to zero when evaluating the viscosity to avoid the anomalies.

Takahashi’s high-pressure correction [30] is used to evaluate binary diffusion coefficients. Since only binary mixtures are considered in this work, the binary diffusion coefficient is the only property needed for the species equations and the evaluation of diffusion coefficients is exact. For more complex mixtures, a mixture-averaged or multi-component treatment of the diffusion coefficient can be adopted [31].

### 2.4 Numerical Schemes

A finite volume approach is utilized for the discretization of the system of equations, \crefeqn:governingEqn:

(8) |

where is the vector of conserved variables, is the face-normal Euler flux vector, is the face-normal viscous flux vector which corresponds to the right-hand side (RHS) of \crefeqn:governingEqn, is the volume of the control volume, and is the face area. A strong stability preserving 3rd-order Runge-Kutta (SSP-RK3) scheme [32] is used for time integration. A Strang-splitting scheme [33] is applied in this study to separate the convection operator from the remaining operators of the system. For LES calculations, the governing equations are Favre-filtered and a Vreman subgrid-scale model [34] is used as turbulence closure.

The convective flux is discretized using a sensor-based hybrid scheme [35] in which a high-order, non-dissipative scheme is combined with a low-order, dissipative scheme to minimize the numerical dissipation introduced. A central scheme which is fourth-order on uniform meshes is used along with a second-order ENO scheme for the hybrid scheme. A density sensor [9, 16] is adopted in this study. Due to the large density gradients present under transcritical conditions, an entropy-stable flux correction technique [16] is used to ensure the physical realizability of the numerical solutions and to dampen non-linear instabilities in the numerical schemes.

Due to the strong non-linearity inherent in the real-fluid state equations, spurious pressure oscillations will be generated when a fully conservative scheme is used [14, 16], and severe oscillations in the pressure field could possibly lead the solver to diverge. A double-flux method [9, 16, 36] is extended to the transcritical regime to eliminate the spurious pressure oscillations. The effective specific heat ratio based on the speed of sound is frozen both spatially and temporally for a given cell when the fluxes at the faces are evaluated. The conservation error in total energy was shown to converge to zero with increasing resolution [16, 36].

## 3 Results and Discussion

In this study, the Spray A configuration [37] is considered, representing a benchmark target of the Engine Combustion Network. The single hole diesel injection is operated with pure n-dodecane at a rail pressure of 150 MPa. Liquid n-dodecane fuel is injected at 363 K through a nozzle with a diameter of 0.09 mm into a 900 K ambient environment at a pressure of 6.0 MPa. The non-reacting case is considered and the ambient gas consists of pure nitrogen. At these conditions, the liquid n-dodecane undergoes a transcritical injection process as illustrated in \creffig:subsupercritical, where the liquid fuel is heated and mixes with the ambient gaseous environment.

A cylindrical computational domain is used in this study with a diameter of 20 mm and a length of 80 mm. The injector geometry is not included in the computation domain. The mesh is clustered in the region near the injector along the shear layers, and stretched in downstream and radial directions. The minimum grid spacing is 4 µm, using approximately 20 grid points across the injector nozzle. The total cell count of the mesh is 4.6 million.

Fuel mass flux and temperature are prescribed at the injector nozzle using the time-dependent rate of injection modeled using the CMT virtual injection rate generator [37], recommended by the ECN with default input parameters for the Spray A case (150 MPa, 6 MPa, 0.0894 mm, 0.90, 713.13 kg/m, and 1.50 ms for injection pressure, back pressure, outlet diameter, discharge coefficient, fuel density and injection time, respectively). Plug flow velocity profile is applied at the nozzle exit without synthetic turbulence. The pressure is prescribed at 6 MPa at the outlet. Adiabatic boundary conditions are applied at all walls. The numerical simulation is initialized with pure nitrogen at 900 K with zero velocity and a pressure of 6 MPa. A CFL number of 1.0 is used during the simulation and a typical time step is about 0.6 ns. The injection process is simulated over a duration of 1.2 ms.

fig:comparison shows a temporal sequence of the injection process from the experiment and current numerical simulation. From this figure it can be seen that the dense fuel jet remains stable in the vicinity of the injector nozzle, for mm. The jet starts to break up downstream of this region and instabilities at the shear layer can be clearly seen. Further downstream, the fuel mixes with the ambient gas and large eddies are present. The spreading angle of the injected fuel is narrow in the injector near-region and becomes wider further downstream where turbulent mixing is present. This is similar to numerical results reported in other studies [8, 10, 11].

The simulation results are qualitatively in good agreement with the experimental images taken by the diffused back-illumination (DBI) method [38, 39] in terms of fuel penetration and overall spreading rate. However, as can be seen from \creffig:comparison, the fuel jet behaves differently in the near injector region. The dense liquid regions are indicated by the blue curves for both the experiment and LES in \creffig:comparison. The liquid region in the experiment is approximately determined by the darkness of the grayscale intensity [39]. A value of 0.6 for the fuel mass fraction is used for the determination of the liquid region in the numerical simulation. It can be clearly seen that the experimental results have a wider spreading angle near the injector. Several possible reasons can explain this discrepancy between experiments and simulation results. One is that the flow inside the injector pipe is not considered in the current simulation. By considering a turbulent inflow profile, stronger shear layer instabilities are expected. Another reason could be due to the limitation of the diffused interface method so that the flow physics in the near injector region is not modeled correctly since interfaces and droplets are completely excluded. However, there is no direct way of determining the liquid region with a diffused interface method. Note that a methodology is developed by Manin et al. [39] to use DBI images to experimentally determine the liquid penetration length, and this methodology is shown to give longer penetration compared to Mie scattering.

The liquid and vapor penetration lengths are extracted from the simulation results using a threshold value of 0.6 and 0.01 for the fuel mass fraction, respectively. The results up to 1 ms after injection are shown in \creffig:pl. The experimental vapor penetration length determined from Schlieren imaging and liquid penetration length from Mie scattering [37, 40] are also shown for comparison. It can be seen that for the vapor penetration an excellent agreement with measurements is obtained. It was found that the utilization of the inflow boundary condition with a time-dependent fuel mass flux is essential for the accurate prediction. Another simulation with a constant mass flux without the initial ramp-up yielded longer penetration length which is consistent with other studies [37]. The sensitivity of the criterion for the determination of the vapor penetration was studied, confirming that the so computed vapor penetration length was insensitive to threshold values between 0.005 and 0.1.

For the validation of the liquid penetration length, a threshold value of 0.6 is used for the fuel mass fraction, providing excellent agreement with the measurements (see \creffig:pl). However, this threshold value is somewhat arbitrary. The sensitivity of the liquid threshold value was also tested, and it was found that for mass fraction values from 0.95 to 0.4 the predicted length varies between 5 mm to 15 mm. It was pointed out in the experimental investigations [38, 39] that the experimentally measured liquid penetration length is sensitive to the choice of the measurement method, the optical setup within one method, the criteria of liquid determination, and the actual geometry of the injector. Indeed, the resolution in either the simulation or the experiment is still not close to fully resolve the interfacial flows near the injector nozzle, if multi-phase flows do exist under these conditions.

The flow structures and mixing behaviors of the injection process further downstream are compared to the measurements of mixture fraction by Rayleigh scattering [41]. The results are shown in \creffig:axial and \creffig:radial. Multiple injections in the experiments provide ensemble-averaged statistics. In the simulation, the statistics of the steady period of injection are obtained by temporally averaging between 0.6 ms and 1.2 ms after the injection. \Creffig:axial shows the comparison of mean and root-mean-square (rms) values of mixture fraction at the center-line of the domain. Good agreement is obtained for both the mean and rms values and the simulation results fall within experimental uncertainties. \Creffig:radial shows a comparison of the radial mixture fraction distribution at three different axial locations (, 25, and 35 mm). As can be seen in \creffig:radial, there is a good agreement in the mean values of the mixture fraction at all three locations, while the simulation predicts slightly higher rms values compared to the experimental data. These results along with the excellent agreement of the vapor penetration length as presented in \creffig:pl, show that the current numerical method is capable of predicting the turbulent mixing process between fuel and surrounding environment downstream of the injector after the dense liquid fuel is fully disintegrated.

## 4 Summary and Conclusions

The need for improved engine efficiencies has motivated the development of high-pressure combustion systems, in which operating conditions achieve and exceed critical conditions. Associated with the transcritical conditions are large thermodynamic gradients as the fluid undergoes mixing and possibly phase transitions. Accurately simulating these real-fluid environments remains a challenge. In this study, a diffused interface method is presented for the modeling of the fuel injection process under conditions relevant to diesel engines. Compressible multi-species conservation equations are solved in conjunction with the Peng-Robinson state equation and real-fluid transport properties. A finite-volume method with an entropy-stable scheme is utilized to robustly and accurately simulate real-fluid flows with strong nonlinearities. A recently developed double-flux model is employed to eliminate spurious pressure oscillations associated with transcritical flows.

LES-calculations are performed to simulate the ECN Spray A target case [37]. Simulation results are analyzed and compared to experimental measurements. The simulation results are qualitatively in good agreement with the experimental images. Different flow dynamics are found in the vicinity of the nozzle where the simulation predicts a narrower spreading angle of the fuel jet. Excellent agreement in vapor penetration length is obtained. The liquid penetration length is found to be sensitive to the selected criterion. Comparisons of the mixture fraction in the turbulent mixing portion of the jet show that the developed modeling capability is able to accurately predict the mixing behavior after the dense liquid jet fully disintegrates.

Observed differences in the flow-field behavior near the injector requires further investigations both numerically and experimentally. Simulations utilizing sharp interface method may provide insights into this question, which is the focus of ongoing work.

## Acknowledgments

Financial support from the U. S. Army Research Laboratory through the Cooperative Agreement W911NF-16-2-0170 is gratefully acknowledged.

### Footnotes

- Corresponding Author: mihme@stanford.edu

### References

- V. Yang, Modeling of supercritical vaporization, mixing, and combustion processes in liquid-fueled propulsion systems, Proc. Combust. Inst. 28 (1) (2000) 925–942.
- W. Mayer, A. Schik, M. Schaffler, H. Tamura, Injection and mixing processes in high-pressure liquid oxygen/gaseous hydrogen rocket combustors, J. Propul. Power 16 (5) (2000) 823–828.
- M. Oschwald, J. J. Smith, R. Branam, J. Hussong, A. Schik, B. Chehroudi, D. Talley, Injection of fluids into supercritical environments, Combust. Sci. Technol. 178 (1-3) (2006) 49–100.
- B. Chehroudi, Recent experimental efforts on high-pressure supercritical injection for liquid rockets and their implications, Int. J. Aerosp. Eng. 2012 (2012) 1–31.
- R. N. Dahms, J. C. Oefelein, On the transition between two-phase and single-phase interface dynamics in multicomponent fluids at supercritical pressures, Phys. Fluids. 25 (9) (2013) 092103.
- R. N. Dahms, Understanding the breakdown of classic two-phase theory and spray atomization at engine-relevant conditions, Phys. Fluids. 28 (4) (2016) 042108.
- D. T. Banuti, M. Raju, P. C. Ma, M. Ihme, J.-P. Hickey, Seven questions about supercritical fluids-towards a new fluid state diagram, in: 55th AIAA Aerospace Sciences Meeting, no. 2017-1106, Grapevine, Texas, 2017.
- G. Lacaze, A. Misdariis, A. Ruiz, J. C. Oefelein, Analysis of high-pressure Diesel fuel injection processes using LES with real-fluid thermodynamics and transport, Proc. Combust. Inst. 35 (2) (2015) 1603–1611.
- P. C. Ma, L. Bravo, M. Ihme, Supercritical and transcritical real-fluid mixing in diesel engine applications, Proceedings of the Summer Program, Center for Turbulence Research, Stanford University (2014) 99–108.
- E. Knudsen, E. M. Doran, V. Mittal, J. Meng, W. Spurlock, Compressible eulerian needle-to-target large eddy simulations of a diesel fuel injector, Proc. Combust. Inst. 36 (2) (2017) 2459–2466.
- J. Matheis, S. Hickel, Multi-component vapor-liquid equilibrium model for LES and application to ECN Spray A, Proceedings of the Summer Program, Center for Turbulence Research, Stanford University (2016) 25–34.
- L. Selle, T. Schmitt, Large-eddy simulation of single-species flows under supercritical thermodynamic conditions, Combust. Sci. Technol. 182 (4-6) (2010) 392–404.
- A. M. Ruiz, Unsteady numerical simulations of transcritical turbulent combustion in liquid rocket engines, Ph.D. thesis, Institut Nationale Polytechnique de Toulouse (2012).
- H. Terashima, M. Koshi, Approach for simulating gas–liquid-like flows under supercritical pressures using a high-order central differencing scheme, J. Comput. Phys. 231 (20) (2012) 6907–6923.
- J.-P. Hickey, P. C. Ma, M. Ihme, S. Thakur, Large eddy simulation of shear coaxial rocket injector: Real fluid effects, in: 49th AIAA/ASME/SAE/ASEE Joint Propulsion Conference, no. 2013-4071, San Jose, CA, 2013.
- P. C. Ma, Y. Lv, M. Ihme, An entropy-stable hybrid scheme for simulations of transcritical real-fluid flows, J. Comput. Phys. 340 (2017) 330–357.
- R. Abgrall, S. Karni, Computations of compressible multifluids, J. Comput. Phys. 169 (2) (2001) 594–623.
- T. Schmitt, L. Selle, A. Ruiz, B. Cuenot, Large-eddy simulation of supercritical-pressure round jets, AIAA J. 48 (9) (2010) 2133–2144.
- S. Kawai, H. Terashima, H. Negishi, A robust and accurate numerical method for transcritical turbulent flows at supercritical pressure with an arbitrary equation of state, J. Comput. Phys. 300 (2015) 116–135.
- R. S. Miller, K. G. Harstad, J. Bellan, Direct numerical simulations of supercritical fluid mixing layers applied to heptane–nitrogen, J. Fluid Mech. 436 (2001) 1–39.
- D.-Y. Peng, D. B. Robinson, A new two-constant equation of state, Ind. Eng. Chem. Res. 15 (1) (1976) 59–64.
- B. E. Poling, J. M. Prausnitz, J. P. O’Connell, The Properties of Gases and Liquids, McGraw-Hill New York, 2001.
- J. F. Ely, H. J. M. Hanley, Prediction of transport properties. 1. Viscosity of fluids and mixtures, Ind. Eng. Chem. Res. 20 (4) (1981) 323–332.
- J. F. Ely, H. J. M. Hanley, Prediction of transport properties. 2. Thermal conductivity of pure fluids and mixtures, Ind. Eng. Chem. Res. 22 (1) (1983) 90–97.
- K. G. Harstad, R. S. Miller, J. Bellan, Efficient high-pressure state equations, AIChE J. 43 (6) (1997) 1605–1610.
- P. C. Ma, D. T. Banuti, J.-P. Hickey, M. Ihme, Numerical framework for transcritical real-fluid reacting flow simulations using the flamelet progress variable approach, in: 55th AIAA Aerospace Sciences Meeting, no. 2017-0143, Grapevine, Texas, 2017.
- T. H. Chung, L. L. Lee, K. E. Starling, Applications of kinetic gas theories and multiparameter correlation for prediction of dilute gas viscosity and thermal conductivity, Ind. Eng. Chem. Res. 23 (1) (1984) 8–13.
- T. H. Chung, M. Ajlan, L. L. Lee, K. E. Starling, Generalized multiparameter correlation for nonpolar and polar fluid transport properties, Ind. Eng. Chem. Res. 27 (4) (1988) 671–679.
- A. M. Ruiz, G. Lacaze, J. C. Oefelein, R. Mari, B. Cuenot, L. Selle, T. Poinsot, Numerical benchmark for high-Reynolds-number supercritical flows with large density gradients, AIAA J. 54 (6) (2016) 1445–1460.
- S. Takahashi, Preparation of a generalized chart for the diffusion coefficients of gases at high pressures, J. Chem. Eng. Jpn. 7 (6) (1975) 417–420.
- R. J. Kee, M. E. Coltrin, P. Glarborg, Chemically Reacting Flow: Theory and Practice, John Wiley & Sons, 2005.
- S. Gottlieb, C.-W. Shu, E. Tadmor, Strong stability-preserving high-order time discretization methods, SIAM Rev. 43 (1) (2001) 89–112.
- G. Strang, On the construction and comparison of difference schemes, SIAM J. Num. Anal. 5 (3) (1968) 506–517.
- A. W. Vreman, An eddy-viscosity subgrid-scale model for turbulent shear flow: Algebraic theory and applications, Phys. Fluids. 16 (10) (2004) 3670–3681.
- Y. Khalighi, J. W. Nichols, S. Lele, F. Ham, P. Moin, Unstructured large eddy simulation for prediction of noise issued from turbulent jets in various configurations, AIAA Paper 2011-2886.
- P. C. Ma, Y. Lv, M. Ihme, Numerical methods to prevent pressure oscillations in transcritical flows, Annual Research Brief, Center for Turbulence Research, Stanford University (2016) 223–234.
- L. M. Pickett, G. Bruneaux, Engine combustion network, Combustion Research Facility, Sandia National Laboratories, Livermore, CA. (http://www.sandia.gov/ECN).
- L. M. Pickett, C. L. Genzale, J. Manin, L. M. Malbec, L. Hermant, Measurement uncertainty of liquid penetration in evaporating diesel sprays, ILASS2011.
- J. Manin, M. Bardi, L. M. Pickett, Evaluation of the liquid length via diffused back-illumination imaging in vaporizing diesel sprays, in: COMODIA2012, 2012.
- L. M. Pickett, S. Kook, T. C. Williams, Visualization of diesel spray penetration, cool-flame, ignition, high-temperature combustion, and soot formation using high-speed imaging, SAE Int. J. Engines 2 (1) (2009) 439–459.
- L. M. Pickett, J. Manin, C. L. Genzale, D. L. Siebers, M. P. B. Musculus, C. A. Idicheria, Relationship between diesel fuel spray vapor penetration/dispersion and local fuel mixture fraction, SAE Int. J. Engines 4 (1) (2011) 764–799.