Analysis of cycle-to-cycle variations of the mixing process in a direct injection spark ignition engine using scale-resolving simulations.
Since the mechanisms leading to cyclic combustion variabilities in direct injection gasoline engines are still poorly understood, advanced computational studies are necessary to be able to predict, analyze and optimize the complete engine process from aerodynamics to mixing, ignition, combustion and heat transfer. In this work the Scale-Adaptive Simulation (SAS) turbulence model is used in combination with a parameterized lagrangian spray model for the purpose of predicting transient in-cylinder cold flow, injection and mixture formation in a gasoline engine. An existing CFD model based on FLUENT vl5.0  has been extended with a spray description using the FLUENT Discrete Phase Model (DPM). This article will first discuss the validation of the in-cylinder cold flow model using experimental data measured within an optically accessible engine by High Speed Particle Image Velocimetry (HS-PIV). Afterwards, the parameterized spray model is validated using experimental data measured in a pressure spray chamber. Finally, results obtained with the combined model are discussed and used to analyze transient mixture formation and to give a detailed insight into cycle-to-cycle fluctuations associated with the turbulent mixing process in an internal combustion engine (ICE).
CITATION: Theile, M., Hassel, E., Thevenin, D., Buchholz, B. et al., "Analysis of Cycle-to-Cycle Variations of the Mixing Process in a Direct Injection Spark Ignition Engine Using Scale-Resolving Simulations," SAE Int. J. Engines 9(4):2016, doi:10.4271/2016-01-9048.
Cyclic combustion variabilities (CCV) are one of the major issues encountered when optimizing gasoline engine combustion. CCV result in high pressure and velocity fluctuations in the combustion chamber, particularly during expansion phase. Figure 1 shows the behavior of the incylinder pressure in an SI Engine and will be used to describe the effects of CCV exemplary. Within a certain operating point with low EGR (black curve) one can see fluctuations of the incylinder pressure. The cycle averaged maximum peak pressure [P.sub.max,mean] is about 41.5 bar, together with a standard deviation [[sigma].sub.p] of about 2.5 bar, which results in a coefficient of variance [COV.sub.pmax] of about 6%. Changing the operating condition to higher EGR, [p.sub.max,mean] stays relatively constant with 40.1 bar, but [[sigma].sub.p] is increasing to 4.1 bar, which results in [COV.sub.pmax] = 10,2%. The same behavior can be seen by looking at the IMEP. The cycle averaged IMEP is almost constant with [IMEP.sub.lowEGR] =6.19 bar and [IMEP.sub.highEGR] = 6.15 bar. Since BMEP was set fixed, differences only occur due to different gas exchange efficiencies. The standard deviation of the IMEP is changing from [[sigma].sub.lowEGR] =0.041 bar to [[sigma].sub.highEGR] = 0.135 bar, while the [COV.sub.IMEP] is increasing from 0.66% to 2.19 %. This does not result in unstable operating conditions for this engine, but it shows clearly the increasing amount of cyclic variability due to different incylinder conditions respectively combustion methods.
Heywood  sets the limit of the [COV.sub.IMEp] to 10%, before serious drivability problems occur. Different groups found unstable operating modes already at lower values of [COV.sub.IMEp] (e.g. 2%  or 4% ).
The origin of this characteristic behavior of CCVs in general goes back to the turbulent in-cylinder flow, developing long before the ignition time, mainly at the intake stage. The different flow fields, developing throughout the intake phase, affect the mixing process and result therefor in cycle dependent initial conditions for the ignition and combustion process. Since fuel efficiency and pollutant emissions both need to be optimized to fulfill new emission standards (EPA Tier 2, CARB LEV-II or Euro 6) with moderate exhaust aftertreatment effort, modern analysis tools must be developed for industrial applications. The coefficient of variance (COV) is an important optimization parameter when trying to fulfill consumer requirements. It can be obtained by CFD simulations when using scale-resolving turbulence models (SRS, in particular the Scale-Adaptive Simulation, or SAS, turbulence model) in combination with already established submodels (e.g., Lagrangian spray models). However, a careful analysis is needed to make sure that all models work properly in combination with each other, since most of the submodels were developed in a classical (U)RANS context.
In this project, a 3D CFD model has been developed to describe and analyze the impact of the turbulent in-cylinder flow on the cyclic fluctuations found during mixture formation in a direct injection spark ignition (DISI) engine. Variations observed concerning global mixture behavior as well as correlations between tumble flow and local fuel/air ratios are then analyzed. Local variations and cyclic fluctuations of the laminar flame speed in the vicinity of the spark plug will be discussed. This is of major importance, since the laminar flame speed controls flame kernel development in the early stage of ignition and combustion . Finally, different regions of interest have to be analyzed carefully and modeled with high accuracy before obtaining a valid CFD model, able to simulate and predict correctly several complete in-cylinder flow cycles including injection and mixture formation.
A key feature of this CFD model is the usage of a hybrid turbulence model, the so called Scale-Adaptive Simulation (SAS) model. Classical turbulence models based on the Reynolds decomposition (RANS) with an eddy viscosity model lead to a calculated turbulent viscosity up to a hundred times higher than that suitable for resolving turbulent structures important for detailed in-cylinder flow analyses . In the past, only cycle-averaged or integral values were considered, for which time or ensemble-averaged simulation results are sufficient. More recently, simulations relying on Large Eddy Simulations (LES; see e.g., [7, 8, 9]) have been increasingly considered in the academic research community. A good and still up-to-date overview is given by . With LES, when using an eddy viscosity based subgrid-scale model, resulting values of the turbulent viscosity are at levels low enough for resolving turbulent structures much smaller than large-scale in-cylinder vortices (responsible for instance for tumble or swirl motion).
However, standard LES simulations are often prohibitive in terms of computing time due to restrictions regarding the numerical mesh and numerical settings, which can influence the stability behavior of the solving process. A number of different alternative scale resolving turbulence models were developed in the past and used for incylinder flow modeling. This includes hybrid models like DES, DDES, IDDES, SDES or SBES, which all use a more or less complex blending function to switch between pure (U)RANS formulation in the near wall zone and LES turbulence modeling in the free stream part of the computational domain based on local grid information and wall distance. These models were successfully used (e.g. by [11, 12, 13]) to capture cyclic fluctuations in the flow field and to show the influence on turbulent combustion.
The need for easy to use turbulence models, which resolve as much turbulence as possible with the given numerical mesh resolution and time step size forced the development of additional global scale resolving turbulence models. Often, these models are not classified as pure LES models, since there is not a direct subgrid scale filter implemented, but the results show, that they are able to resolve nearly as much turbulence, compared to pure LES models. [14, 15] used a so called Partially-Averaged Navier-Stokes Simulation (PANS) model, which is originally based on a RANS k-[epsilon]-model. The eddy viscosity is calculated using additional transport equations for unresolved turbulence quantities (kinetic energy and dissipation rate). This approach was extended by  and was converted to the transport equation of the RANS k-[zeta]-f model. The adjustment of the eddy viscosity to resolve a higher amount of turbulence is done by using the ratio of grid size and integral length scale. This results in a flow field dependent value for the eddy viscosity.  showed, that PANS is able to produce more accurate flow fields and swirl numbers in comparison with a RANS-k-[zeta]-f model in a portflow application. Simulation results were compared with PIV and flow bench measurements.  showed, that the PANS model is able to resolve a certain part of the turbulence, but show an instable behavior in an application in real engine geometries. This behavior was explained by the highly decreased resolution parameter [f.sub.k], which leads to very small eddy viscosities in the valve gap.  developed a new scale resolving turbulence model, the Dynamic Length Scale Resolution Model (DLRM), which dynamically adapts the eddy viscosity as a result of a comparison between modelled and potentially resolvable length and timescales. The latter are estimated using grid and time step size information. In  a simplified engine geometry with a fixed, central positioned intake valve was used as a validation test case. DLRM simulations on a coarse (700 K cells) and fine (4 M cells) grid size were compared to Large Eddy Simulations using a WALE subgrid scale model on a 12 M cells mesh and with LDA measurements. The DLRM results showed a remarkable agreement in averaged velocities and rms fluctuations, partially with better agreement with the experimental data than the LES-WALE model.
An interesting alternative for predicting the transient in-cylinder flow in a similar context is the Scale-Adaptive Simulation (SAS) turbulence model. This turbulence model was first described by Menter [18,19, 20] in 2003. It was shown that SAS delivers additional, fine-scale flow information like for a LES turbulence model, and leads to improved comparisons with experimental results compared to classical RANS models, like k-[omega]-SST. Previous studies with geometries close to real ICE confirmed the promising outputs obtained with SAS [12,21].
The given examples regarding the usage of dynamically adaptive scale resolving turbulence models show, that they are interesting and already established alternatives to classical LES models. Time depended meshes and the time step sizes used in incylinder flow simulations are not always sufficient for pure LES models. The adaptive behavior of the discussed adaptive SRS models, regarding given mesh and time step sizes, makes them favorable for scale resolving incylinder flow simulations.
The final goal in this study is to develop a CFD model and a modeling and analyzation strategy that allows detailed analyzation of cyclic combustion variability for industrial applications. This includes the usage of an adaptive SRS model, which will be validated against PIV measurements. This ensures the prediction capabilities of the chosen SRS model in the specific case of this study.
Scale-Adaptive Turbulence Model
Most RANS-based turbulence models use an eddy viscosity approach, where the influence of the non-resolved turbulent motion on the time-or ensemble-averaged flow field is described using an assumption of turbulence homogeneity and isotropy. The unknown Reynolds stress tensor resulting from the Reynolds decomposition is then calculated based on the (known) mean velocity gradient tensor and on the eddy viscosity:
[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (1)
A common approach to estimate the eddy viscosity is to use at least one variable describing turbulence, which is calculated via its own transport equation. Nevertheless, most turbulence models currently employed for engine simulation use two-equation models, namely k-[epsilon], k-[omega]-SST or k-[zeta]-f. The SAS turbulence model used in this study is originally based on the k-[omega]-SST model , so that the eddy viscosity is calculated following:
[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (2)
Here, [alpha]* is a low-Reynolds number correction and results in unity for high-Reynolds number flows; [[alpha].sub.1] is a model constant ([[alpha].sub.1]=0.31), [F.sub.2] is a blending function of the SST extension and S is the magnitude of the strain rate.
The difference between the standard k-[omega]-SST and the SAS model is a source term in the transport equation for the turbulence dissipation frequency ([omega]-equation). This important term involves the ratio between resolved and unresolved turbulent length scale, based on the von Karman length scale. The latter is derived from an exact transport equation describing the turbulence length scale. Further details are described in [18,19]. Finally, the source term is defined as follows:
[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (3)
with the model constants [[ETA].SUB.2]=3.51, [[sigma].sub.[PHI]]=2/3 and C=2. The integral length scale L is defined as
[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (4)
with the model constant [c.sub.[micro]] = 0.09. The von Karman length scale [L.sub.vK] is calculated via
[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (5)
The comparison of both length scales in this model results in a reduction of the eddy viscosity, when the formation of a resolved turbulent spectrum is detected by the SAS model, for small values of the von Karman length scale. Therefore, the model can act in a "scale-resolving" or "pseudo-LES" mode. This is particularly true for strongly unsteady flows, where fast variations in time can be introduced via boundary or initial conditions, or can develop spontaneously through non-linearity effects appearing in the Navier-Stokes equations in combination with complex geometric features. In the case of ICE, flow separation at the intake valves can trigger the model . The potential of SAS for in-cylinder flow simulations has been already confirmed by different research groups [12,21].
Engine Model - Cold Flow
The simulation domain is a full four-valve pent-roof combustion chamber with the characteristic bore of a passenger car. It incorporates the intake and exhaust ducts, which are numerically resolved up to 180 mm upstream the intake valves and 80 mm downstream the exhaust valves. Additional information is given in table 1.
The size of the numerical mesh employed in this study (see Figure 2) is 3.6 million cells in TDC and 6.3 million cells in BDC. The ICE Tool of the ANSYS Workbench was used for meshing after adapting it to the problem at hand. A hybrid meshing approach was used for the dynamic mesh in the CFD solver. To improve the accuracy of the prediction concerning flow separation at the intake valve, the mesh density was increased in this region. The average cell size was set there to 0.2 mm. In order to resolve a sufficient part of the turbulent spectrum, the cell size in the combustion chamber was set to 0.6 mm. The maximum cell size was limited to 2 mm. This value was only reached in the intake and exhaust ducts.
The boundary conditions for the CFD simulation were taken from experimental data of individual cycles measured in the intake- and exhaust ducts of the test-bench engine, 120 mm upstream the intake valves and 120 mm downstream the exhaust valves (see figure 3). The spatial distances between these measurement locations and numerical boundaries were found to have only minor influence on the incylinder flow field. Assuming that the pressure wave moves with the speed of sound (roughly 340 m/s), the spatial distance between the numerical boundary and the real location of the pressure signal measurements results in a shift of the pressure boundary signal of 1.5 CAD at 2000 rpm. The benefit of the chosen position of the inlet boundary is an increased duct volume, so that the influence of the scavenging air and resulting EGR (in combustion simulations), based on camshaft timings, can be directly resolved.
For each simulated cycle, specific boundary data have been used in terms of pressure curves (see figure 3). The COV of the intake and exhaust pressure curves are always below 0.5 % (see table 2). In combination with SAS, this modeling strategy leads to cycle-dependent flow fields as simulation results. Preliminary (U)RANS calculations with a k-[omega]-SST turbulence model showed only minor differences in the cylinder flow field. It can be stated, that the differences in the cycle dependent boundary conditions are not high enough to overcome the effects of the high turbulent diffusion of a common (U)RANS model. Since no turbulent fluctuations are introduced at the inlet boundary, an additional benefit of these kind of boundary conditions is simple implementation, when measurements are available.
Temperatures for inlet and outlet boundary conditions were measured 180 mm upstream the intake valve and 180 mm downstream the exhaust valve. The measurements were not crank angle resolved, so that the inlet and exhaust boundary condition were set as constant values. Surface temperatures at the different boundaries (e.g. intake and exhaust ducts, cylinder head, piston, liner) were set constant with values between 303 K and 383 K. Initialization was done by first simulating a preliminary cycle with the use of a RANS-k-[omega]-SST turbulence model. Starting from this point, a different amount of multicycle simulations were started, each consisting of six consecutive cycles, whereas the first cycle was always dropped to ensure independence from initially prescribed values. The amount of the multicycle simulations is different for the analyzed operating points, resulting in different amounts of total cycles to analyze (see Table 2). For operating point OP 2, a total number of 30 cycles was reached by using six multicycle simulations. The dynamic mesh approach of FLUENT was used to describe the moving boundaries (piston and valves).
Multiple cycles of three different non-fired operating points have been simulated, compared and validated with experimental data measured at the test-bench engine. Data for the validation process include incylinder pressure data measured in the combustion chamber, and cycle-averaged as well as rms values of the velocity, measured with High-Speed Particle Image Velocimetry (HS-PIV). Simulation and experimental results from HS-PIV measurements have been interpolated on congruent 2D post-processing grids for accurate analysis.
Three different operating points have been simulated with several individual cycles in order to gather statistics, focusing mainly on operating point 2. More details can be found in Table 2. The computational time needed for one cycle can be estimated by 4600 CPUh on Intel Xeon-E5 Cores. Comparisons with the same setup using DDES and WMLES turbulence models showed, that the convergence performance of the SAS model overcompensated the additional computation time for the additional transport equations respective more complex turbulence modeling procedure. Since convergence behavior and simulation stability was a major aspect, the decision was in favor of the SAS model.
Spray Model for Fuel Injection
For simulating injection of gasoline fuel the FLUENT Discrete Phase Model (DPM) was used. This is a Lagrangian approach classically employed for modeling of diluted liquid or solid secondary phases. For a detailed description, see [23, 24, 25].
To model injection, spray breakup, and evaporation, different submodels were used. The injection process is modeled by using a blob injection model with a uniform diameter equal to the nozzle hole diameter. Primary and secondary breakup were modeled using the KHRT model . This model considers Kelvin-Helmholtz instabilities at the surface of the liquid core of the spray to initialize breakup. This results in the generation of small droplets and ligaments forming out of the initial blobs. Additionally, Rayleigh-Taylor instabilities are modeled, generating smaller droplets in the diluted part of the spray. Turbulence production or damping at subgrid scale through droplets smaller than the grid size is considered using an implemented model described in .
In all simulations, a single-component surrogate was used instead of a multi-component fuel, which is a major but often used simplification. Gasoline is a mixture of more than hundred single hydrocarbons, which all have different vapor pressures. To describe evaporation , iso-octane was used as single component, since its vapor pressure lies near the average vapor pressure of the main components of gasoline fuel (see Figure 4). The properties of the fluid are computed using standard databases .
The injected mass flow-rate was evaluated experimentally. The injection velocities were calculated using mass conservation in the injector nozzle and assuming an incompressible fluid. The validation of the spray model was done using experimental data from measurements at an optically accessible pressure chamber. The geometrical features of the chamber were used to setup up a CFD model for the gasoline injection process. Boundary conditions for temperature and initial values for the chamber pressure and temperature were set based on measured experimental values. Since there is no data available from the experiment regarding initial turbulence values in the chamber, values from a similar case  were applied. The numerical results were then compared with measurements. Available experimental data include High Speed Photography (HSP) and Phase-Doppler Anemometry (PDA) measurements delivering particle size distributions at 40 mm distance from the nozzle hole. The HSP data was post-processed to obtain the tip penetration length over time. Different operating points were analyzed, representing a wide range of possible engine operating points. In order to predict spray penetration and evaporation correctly, three operating points were used to fit the model, as described in Table 3.
Validation of the Engine Cold Flow Model
This first validation was done by comparing simulation results and measurements of the in-cylinder pressure and the velocity field in different planes within the combustion chamber.
As seen in Figure 5, there is a good agreement between prediction and measured in-cylinder pressure. Both results match very well during the gas exchange and early compression phase. In the timespan near the top dead center (TDC), there is a slight difference in comparison with the experiment. This can be explained by a loss of mass (blow-by) in the experimentally employed combustion chamber. At higher engine speeds, blow-by decreases. The presented in-cylinder model is optimized for operating point 2 (2000 rpm / 1 bar intake pressure), by adjusting compression ratio through a shift of the piston position to obtain the same peak pressure in TDC. The slight overshoot of pressure in the simulation results at TDC in OP 1 and the corresponding undershoot in OP 3 can be explained by the different engine speeds, since blow-by is a function of time, not of crank angle. The blowby of the real engine is decreasing with higher engine speeds, so that the pressure difference in comparison with the simulation is negative for high engine speeds, and vice versa for lower engine speeds. Blow-by is, in addition to wall heat transfer, also the reason why the maximum pressure is shifter to earlier crank shaft positions in the experiment. However, and in agreement with the scientific literature , the impact of blow-by on the velocity field is assumed to be minimal, so that blow-by is not considered at all in the CFD model.
The effect of blow-by can also be observed during the late expansion phase, when the measured in-cylinder pressure is lower than the volume-averaged in-cylinder pressure in the simulation. After the exhaust valves open (EVO), the pressure in the simulation drops to levels nearly identical to the measurements. In the gas exchange phase the pressure curves match again very well. At IVC the averaged pressure deviation between simulation and experiment is below 0.025 bar, which indicates, that nearly the same mass is contained in the combustion chamber. Possible differences concerning the really trapped mass can occur due to missing information concerning temperature at the boundaries. Since the optical engine has fixed operating temperatures for engine oil and cooling water, wall temperatures values were chosen based on in-house calculations and previous experience. Considering the obtained agreement, it is confirmed that this procedure works very well.
Comparisons of the cycle-averaged velocity fields and of the rms values of velocity were carried out next. In order to understand the experimental procedure, a short description of the setup at the transparent test bench engine (see figure 6) is given first.
The High-Speed PIV System consists of two lasers (Quantronix Hawk II and Quantronix Hawk Pro) with an output wavelength of 532 nm. The images were captured with a LaVision HighSpeedStar 5 camera with a resolution of 512 x 512 pixel. The acquisition frequency was 2640 Hz at OP 1/OP 2 and 3960 Hz at OP 3. The time delay between one double image varied between 5 and 15 [micro]s. Peanut oil was used as seeding particles, which delivers best performance regarding atomization, flow following and vaporization behavior compared to other seeding substances like plant oil, magnesium oxide and propylene carbonate. The smoke number is above 433 K. In combination with a relatively low compression ratio of 7:1, vaporization occurs only at the time span around TDC. The seeding density was still high enough to calculate velocity vectors. The peanut oil was introduced and atomized 1 m upstream of the cylinder head. The raw images were post-processed with the commercial software LaVision Davis 8.16. For validation, images were analyzed every five CAD in 74 consecutive cycles from 450 CAD b. TDC to 30 CAD a. TDC.
The engine setup allows measurements as either front or side view through the glass liner, or as bottom view through the pancake piston. A schematic illustration of the optical setup can be seen in figure 7.
A representation of the two different planes finally used for validation in the present study can be seen in figure 8. By selecting these planes, essential flow features, in particular the tumble motion, are captured.
Comparisons will be shown for each plane at three different crank angle positions, namely 270, 180 and 90 CAD b. TDC for operating point OP 2.
The comparison of the cycle-averaged values of velocity magnitude at different crankshaft angles shows overall a good agreement. Figure 9 depicts the mean velocity field calculated with the two components of the velocity vector tangential to the vertical cut-plane through the middle of the combustion chamber at OP 2 for different crank angle positions. Since near-wall PIV measurements are very difficult due to laser reflections, the corresponding, "low-confidence" area has been cut out in the images. Even in the remaining part of the domain, the margins correspond to a lower signal-to-noise ratio must thus be considered with caution.
In the early intake phase (270 CAD), a slight overprediction of the velocity magnitude is observed in the simulation in the region of the intake jet. This can be found for all analyzed operating points. This phenomenon can be caused by the high turbulent fluctuations of the two incoming air jets through both intake valve gaps. In the experiments, this can lead to seeding particles with very high velocities normal to the laser measurement plane, reducing measurement quality. It is then possible that many fast particles flow through the plane without being tracked. Overall, the comparison between SAS and HSPIV is still very good.
Additionally, it can be seen, that the tumble vortex starts developing in the middle of the observation window. At 180 CAD BTDC the upper side of the large tumble motion is clearly visible in both simulation and experiment. The agreement concerning tumble center is considered as a major quality criterion when comparing simulation results with experimental data. Here, the tumble core is estimated to be in the bottom right part of the cylinder, both for SAS and HS-PIV, but cannot be captured by the PIV, because of the limited window size and signal quality in the lower part of the cylinder. An analysis of the tumble trajectory (see figure 14) shows, that the location of the tumble core is converging into the bottom right corner towards the end of the intake phase for both simulation and experiment.
A comparison of the mean velocities at 90 CAD BTDC shows again that there is a very good match regarding the global flow field and the tumble core position.
In Figure 10 corresponding flow fields are presented in the valve plane. One can clearly see the incoming air jet in the SAS Simulation in the top left corner at 270 CAD BTDC. In the experiments, the useable illumination plane in the valve cross-section is smaller than the middle plane due to the non-centered plane position in the combustion chamber and increased laser reflections. As a consequence, the incoming jet is not clearly visible in the PIV results. Low-quality experimental data have been obtained in the top-left corner at all crank angle positions for the valve plane. Nevertheless a comparison of the flow field in the center of the observation window shows very good agreement. Again, tumble core position and evolution as function of the crank angle positions agree very well.
For further comparison, the local velocity components at eight locations inside the combustion chamber has been compared (see figure 11 and 12)
Figure 11 shows the cycle averaged velocity components tangential to the measurement plane and the rms fluctuations in four different locations at the mid plane. The SAS datasets were calculated using 30 cycles, while the experimental data was obtained using 74 cycles. In location 1, the averaged velocity is over predicted by the simulation, which can be seen by looking at the u- and v-velocities (first row).
This over prediction occurs in the first part of the intake phase until 200 CAD BTDC. The same behavior can be seen at location 2 regarding the v-velocity. All other locations show only minor differences between simulation and experiment, so that the overall matching was found to be good. The characteristic flow behavior can be reproduced. By looking at the bottom eight graphs of figure 11 it can be seen, that the fluctuations of the velocity can be captured very well. There is a slight over prediction of rms fluctuations in location 2 in the early intake phase, which decreases until 180 CAD BTDC to the level of the experiment. The overall level of velocity fluctuations in all other points is very similar. At the end of the cycle, from 60 to 40 CAD BDTC on until the piston reaches the specific probe location, there is an increasing of rms fluctuations visible in the experiment. This can be seen especially well in location 1 and 2 and partially in location 3 and 4. This was found to be the result of the piston motion, which decreases the signal quality, so that higher fluctuations of the velocity components start appearing.
Figure 12 depicts again four locations, only with an offset onto the valve plane (see figure 8). One can see minor differences in the mean u- and v-velocities in location 5 and 6, again with a slight over prediction in the simulation in the early intake phase. The location of the intake jet is not exactly the same, so that the vertical and horizontal velocity components are different. There is a convergence regarding the velocity components visible, which results in a good matching after the early intake phase. In location 7 and 8, the differences in the intake phase are smaller. Again, the rms fluctuations look very similar. In location 5 and 6 the simulation show higher amount of turbulence, which is a result of the higher level of large scale fluctuations in the simulation, which could not be captured in the experiment. At location 7 and 8 the behavior of the turbulence intensity is the same. The peaks in the experiment are again a result of the lower PIV signal quality.
Overall it can be seen that the mean velocities can be predicted well by the simulation, especially in the later intake and compression phase. A high amount of cycle-to-cycle velocity fluctuations are properly resolved by SAS, and show features similar to the experimental data.
A comparison between modelled and resolved turbulent kinetic energy (TKE) shows that in the later intake phase and during compression phase the amount of resolved TKE is above 80% in major parts of the combustion chamber (see figure 13). The presented resolved TKE was calculated using the ensemble averaged velocity fields at fixed crank angle positions as mean velocites. Since it is hard to distinguish between large scale fluctuations and small scale turbulent fluctuations, both of them are included when calculating the velocity fluctuations based on phase-dependent data. Nevertheless, the amount of resolved fluctuations is far beyond values characteristic in (U)RANS simulations.
The upper part of the combustion chamber consists of a tetrahedral mesh, which has a slightly smaller averaged cell volume than the hexahedral mesh in the lower part. The different mesh sizes cause different levels of resolved velocity fluctuations, because of lower turbulent viscosities. Due to the high velocity magnitude, the fluctuations of the intake jet are not resolved that much, as it can be seen at 240 CAD b. TDC. This is also the reason, why the comparison with PIV measurements is important. Due to the capability of the SAS model, to switch smoothly between (U)RANS and scale resolving mode, one aspect of the chosen modeling approach was, to determine the quality of the simulation, when some parts of the turbulent flow are not highly resolved, but modeled more in an (U)RANS context. This can result in major savings of computational resources, since the timestep size can be chosen moderately.
The final conclusion is, that a major part of the turbulence can be directly resolved by SAS simulations in this DISI engine model setup.
Because the mixing process is mainly controlled by the charge motion in the combustion chamber, a comparison of the tumble movement was performed. This main vortex structure is created during the intake phase and moves through the combustion chamber due to piston motion. This can partly be seen in Figures 9 and 10. The vortex center identification method called [GAMMA]1-Criterion  was used to capture this movement based on the velocity field. A comparison of the trajectory of the tumble vortex showed good agreement, demonstrating that this important flow feature can be predicted with the CFD model (see figure 14). Since the tumble vortex is not a perfectly circular flow feature, the areas of high values of the [[GAMMA].sub.1]-Criterion has been displayed in figure 14 in addition to the tumblecore positions, which were calculated using the global maximum of [[GAMMA].sub.1] In the intake phase between 280 and 260 CAD b. TDC the tumble movement is created and stretched, due to the incoming air jet. This is displayed in figure 14 by the non-circular areas of [[GAMMA].sub.1] > 0.5. This is also the timespan, when there is a discrepancy between simulation and experimental data, which was discussed earlier. At 180 and 200 CAD b. TDC, the tumble positions converge to each other in the bottom right corner. Due to the low signal to noise ratio, a clear determination via [[GAMMA].sub.1]-Criterion was not possible between 200 and 140 CAD b. TDC. The window size of the optical accessible engine was limited to 55 mm below the cylinder head gasket. In the compression phase, between 140 and 40 CAD b. TDC, the tumble moves upwards to the exhaust valves. Between 60 and 80 CAD b. TDC the tumble core moves slightly towards the intake valve side before the breakdown begins. This behavior is visible in the experimental data and also in the simulation. The same behavior is described in , who compared his LES results with experimental data. The breakdown of the tumble is of major importance to create turbulent kinetic energy at small scales, supporting mixing, ignition and turbulent combustion in the gasoline engine. Using the [[GAMMA].sub.1]-criterion, an estimation of the breakdown time could be obtained, around 60 CAD BTDC, since the size of the area with [[GAMMA].sub.1] > 0.5 is decreasing. This is in compliance with the characteristic decrease of the tumble ratio. This time was identically found in the simulation and in the experiment.
A similar behavior can be seen in the valve plane (see figure 14 bottom two graphs), where the tumble moves towards the bottom right corner in the intake phase, then moves upwards to the exhaust side in the compression phase, before it finally moves towards the inlet side at the end of compression. Figure 14 shows, that the characteristic movement of the tumble vortex can be predicted within a certain accuracy. It is the author's conclusion, that the cold flow model is sufficiently validated at this point.
Validation of the Injector Spray Model
A comparison of the simulated transient spray with experimental data shows overall a good agreement. Figure 15 presents a comparison of spray shape in experiment and simulation for operating point OP s3 (see Table 3) at different times. The simulated penetration and spray angle is in good agreement with the measured data from High Speed Photography. There is a difference at the beginning at 0.5 ms after start of injection. At this time, the penetration of the spray is higher in the experiment. In addition to aerodynamic breakup processes, flash boiling occurs, which can form an increasing penetration at the beginning of the injection . From the literature and process conditions, flash boiling is only exists in OP s2 and OP s3. This phenomenon happens in situations where fluids experience a sudden pressure drop, so that the aggregate state at thermodynamic equilibrium is changing from liquid to gas in an extremely fast manner. This results in sudden evaporation with additional spray breakup. The critical chamber pressure, at which flash boiling could be observed at the stated temperature ranges, was at 1.2 bar and lower pressures. The KHRT model is not considering this kind of additional break-up effects. To consider indirectly the described effect in our simulations, the KH-RT constants and initial spray angle have been chosen to represent the sum of all spray breakup phenomena (spray penetration, measured droplet size distributions and spray shape) as good as possible. In this study, only a single injection will be used as injection strategy for the engine model, resulting in injection durations of 2 ms and longer. Therefore, the main concern lies within particle size distribution and spray penetration in the main phase to model the evaporation and mixing with accuracy.
The penetration converges towards the experimental data after 1 ms, so that the main injection phase is modelled correctly. The PSDs measured at 40 mm from the nozzle hole are in close agreement in experiments and simulation (see figure 16). This is important, since there is a direct impact of the droplet size distribution on the penetration length. The agreement concerning PSD is thus a major advantage. The initial droplet radius of the blob injection decreases over space and time due to breakup. Getting the correct PSD, the volume-to-surface ratio of the droplets will be correctly reproduced, leading to correct evaporation rates and, ultimately, mixing.
An analysis of the time step size and mesh dependency was done by changing grid size, mesh type and time step size to values, which are used in the engine model (see table 4 and figure 17). In the simulation setups 2, 3 and 4 the changes for penetration length and spray angle are below 3%. Also only minor changes in PSD occur, so that the engine setup was found to be adequate. The results are confirmed by similar investigations . This shows that the main spray properties are well described by the employed numerical model, so that mixing can now be considered.
Mixture Formation in a DISI Engine
Some changes had to be done after combining cold flow and spray models, before predicting the transient mixture formation in the gasoline engine geometry. The main modification is the adjustment of the total injected fuel mass, resulting in a global stoichiometric fuel-air-ratio of unity. The start of injection was set constant to 280 CAD b. TDC in all operating points. Since the injection pressure was also assumed to be constant, the injection rate profile was modified by increasing the duration of the maximum injection rate. This results in different injection durations in terms of crank angle degree in all operating points. Taking the mass fraction of oxygen in air to be 0.2315, the air requirement for iso-octane is calculated to be 15.16 [kg.sub.air]/[kg.sub.fuel] Corresponding boundary conditions of the engine injection case can be taken from Table 5.
In Figure 18 the local distributions of the equivalence ratio are shown, for three individual cycles at operating point [OP.sub.mixture] 2. High cycle-to-cycle fluctuations are observed concerning the fuel distribution near the spark plug, but also in the whole combustion chamber. Additionally, small-scale turbulent fluctuations are visible, which could not be observed in classical RANS simulations. At 80 CAD BTDC the increased mixing process due to tumble vortex breakdown starts. Areas corresponding to exceedingly rich and lean mixtures reduce continuously from 80 CAD BTDC until TDC.
The cycle-to-cycle fluctuations of the equivalence ratio start with the injection process. Since the injection parameters of the spray model are identical for each cycle of the same operating point, only the interaction with the turbulent flow field can lead to a different spray behavior. As a consequence, a different evaporation and mixing process is observed. In the timespan between the end of injection and top dead center the formation of a homogeneous mixture is not completed overall, confirming observations from the literature (e.g. ).
To analyze the mixing behavior of different cycles, probability density functions (PDF) have been computed. The mass-weighted values of the equivalence ratio (ER) in each cell of the combustion chamber were binned. The resulting values were then normalized, to ensure that the integral of the pdf is equal to unity, in order to generate comparable plots. After that, all cycles of the analyzed operating points were sorted by comparing the mixture quality at TDC. Mixture quality was quantified using the height and width of the PDF at TDC for each cycle. To analyze the differences in mixing behavior, the PDF of the two cycles with the best mixture formation were compared with the worst two cycles. Figure 19 shows corresponding results at operating point [OP.sub.mixture] 2 and reveals very large differences. The PDFs of the best cycles (cycle 6 and 15, shown at the top) correspond a very homogeneous system close to stoichiometric conditions. While the behavior at TDC on the lean side (ER between 0.8 and 1) is similar for all cycles shown in Figure 19, the two worst mixing cycles (cycle 9 and 28, presented at the bottom) show considerably differences at high equivalence ratios. Considering an ER of 1.15 as a threshold, the PDFs of cycles 6 and 15 are almost 0 on the rich side of the limit, while cycles 9 and 28 show probabilities between 0 and 0.02 for ER > 1.15. Finally, comparing the peak values near stoichiometry clearly reveals that such conditions dominate much more in cycles 6 and 15 compared to cycles 9 and 28.
Looking at the time span of the tumble vortex breakdown from 80 CAD BTDC to TDC it is clearly visible that turbulent diffusion leads to a fast mixing process. The PDFs at 80 CAD BTDC in Figure 19 reveal for all cycles lean cells in the range of ER between 0.5 and 0.7. However, many rich cells with an ER > 1.5 are only found in the worst mixing cycles, 9 and 28, while they are not observed in the best mixing cycles. Therefore, it appears that vortex breakdown induces strong mixing and is able to compensate for lean conditions, leading to stoichiometry until TDC. At the same time, vortex breakdown cannot sufficiently solve the problem of rich pockets, so that exceedingly rich conditions found at an early stage are still noticeable at TDC.
A different mixing behavior can be seen in operating point [OP.sub.mixture] 1 and [OP.sub.mixture] 3 (see figure 20 and 21), where the amount of injected fuel mass is nearly doubled, compared to [OP.sub.mixture] 2. Due to longer injection timings the evaporation lasts longer during the compression stroke, resulting in pdf locations on the lean side at 80 and 40 CAD b. TDC. The evaporation process accelerates as a result of the increasing gas temperatures in the final compression phase, until it is nearly completed at TDC. Therefor the pdf at TDC is right around the stoichiometric value of unity. By categorizing the cycles similar to [OP.sub.mixture] 2, cycle 5 and 1 are the corresponding best mixing cycles in [OP.sub.mixture] 1 and 3. A comparison with the both worst mixing cycles (cycle 2 in [OP.sub.mxture] 1 and cycle 8 in [OP.sub.mixture] 3) shows major differences in mixing already at 80 CAD b. TDC on the lean side. In [OP.sub.mixture] 1 cycle 2 has a peak at ER = 0.65, while cycle 2 has only a moderate values around that ER. A similar situation can be seen in [OP.sub.mixture] 3 at ER = 0.74. These differences in mixture formation between cycles are higher than the differences observed in [OP.sub.mixture] 2 discussed earlier. The convection and diffusion processes in [OP.sub.mixture] 1 and 3 are not strong enough to overcome the mentioned deficits, resulting in a wider pdf at TDC for cycle 2 ([OP.sub.mixture] 1) and cycle 8 ([OP.sub.mixture] 3).
Showing now the mixing process by looking at isosurfaces of the equivalence ratio at one specific crank angle position, the major differences can be visualized for [OP.sub.mixture] 2 exemplary (see Figure 22). At 80 CAD BTDC, rich zones are mostly found at the outer boundary of the tumble vortex, and are much predominant in the two worst cycles, 9 and 28 (bottom row of Figure 22). This points out the importance of modern injection strategies like split injection or geometry-optimized targeting. With a single injection, there is a high probability that only a part of the combustion chamber is filled with the injected fuel.
To quantify the cycle-to-cycle fluctuations near the spark plug gap, two different regions were analyzed, a small sphere [S.sub.1] with a radius r = 1 mm and a large sphere [S.sub.2] with r = 10 mm (see figure 23), both centered at the spark plug gap. Using [S.sub.1], small-scale fluctuations that directly affect the ignition and early flame kernel development within the spark plug gap can be analyzed. The influence on such fluctuations on flame growth can be estimated using [S.sub.2]. Since literature considers that overall cycle-to-cycle fluctuations mainly originated from the early stage of combustion , this analysis is particularly interesting.
Figure 24 shows the local, mass-weighted, volume-averaged equivalence ratio, mixture temperature and resolved turbulent kinetic energy calculated in the presented two spherical volumes at 20 CAD BTDC. It is clear that the local fuel-air ratios in the vicinity of the spark plug fluctuate strongly from one cycle to another.
It can be seen that the volume-averaged equivalence ratio of [S.sub.1] is varying between 0.82 and 1.13, whereas fluctuations in the sphere [S.sub.2] are smaller, which is of course a direct result of the higher amount of cells used in the volume-averaging process. Nevertheless, the cyclic fluctuations are still noticeable in the larger volume, showing values of the equivalence ratio between 0.91 and 1.04 as function of the cycle.
The variation in the local temperature is a direct effect of the evaporation process. In areas with a high evaporation rate, the internal energy of the air surrounding the liquid droplets is used to overcome the latent heat of iso-octane. Species mixing and heat diffusion is mainly controlled by the turbulent diffusion. Therefore, the correlation between rich (respectively lean) areas and lower (resp. higher) temperatures is very high. This is a direct result of the turbulent transport model, since both transfer processes are described using a gradient-diffusion hypothesis. The effective diffusivity is the sum of the molecular and the turbulent diffusivity. The latter is calculated using the turbulent viscosity, the turbulent Prandtl number [Pr.sub.t] (for heat diffusion) and the turbulent Schmidt number Sc (for species mixing). [Pr.sub.t] and Sc are constant, equal to 0.85 and 0.7, so that diffusion of heat and species take place in a similar manner. Nonetheless, the resulting temperature variations shown in figure 24 have an additional impact on ignition and flame kernel growth.
To give additional information about the resolved velocity fluctuations that correlates with the resolved turbulent kinetic energy (TKE), the normalized quadratic velocity fluctuation (NQVF) was calculated using the mass-weighted volume-averaged velocity components and the mass-weighted values in each cell of the presented spherical volumes. It is clear that NQVF shows very high cyclic fluctuations. This underlines again the ability of SAS to capture flow field fluctuations, thus directly affecting the turbulent combustion process after the flame kernel surface becomes large enough to get wrinkled by turbulent eddies.
In addition to figure 24, the coefficients of variance for the equivalence ratio, temperature and NQVF has been displayed in table 6. The ratios of the resulting COV for the equivalence ratio and temperature of both spheres are nearly identical. This confirms the strong coupling between mixing and temperature diffusion as discussed earlier. The COV for the NQFV is nearly identical for both spheres. The high cyclic fluctuations of the velocity vector in both spheres are direct results of the interaction with the spark plug geometry. The level of [COV.sub.NQVF] shows, that cyclic variabilities are nearly equal in both areas.
As a first step towards quantifying the possible impact of cyclic variations of fuel/air mixture and temperature on the combustion, the changes in laminar flame speed resulting from the fluctuations have been computed. For this purpose, a detailed reaction mechanism was used for iso-octane . The calculation of the laminar flame speed for the given specific equivalence ratios and temperature levels was carried out using CHEMKIN. Results are shown in figure 25. The resulting COV of the laminar flame speed is also displayed in table 6. It is clear, that the ratio of COV between both spheres are at the same level as the ratios of the equivalence ratios and temperatures.
Comparing the range of resulting laminar flame speeds in [OP.sub.mixture] 2 for all 30 cycles, one can conclude that variations of [S.sub.L] between 40 cm/s in the leanest cases and 55 cm/s in richest conditions must be expected. Due to the strong correlation between early kernel growth and laminar flame speed, relative differences in flame kernel growth by about 30 % will occur just after ignition. Since this analysis does not take into account the impact of aerodynamics and turbulence, even larger fluctuations have to be expected.
A numerical model for the simulation of the cold flow, injection and transient mixture formation in a real DISI engine has been developed using a commercial CFD software. Different submodels were presented, discussed and validated by comparison with experimental data. An analyzation strategy was presented, based on the calculation of probability density functions and the analyzation of equivalence ratio, temperature and velocity fluctuations near the spark plug. Calculations of the laminar flame speed at the conditions found in the simulation were combined with simulation results to estimate fluctuations of the early flame kernel development.
The functionality of the SAS turbulence model within a CFD model of a real ICE geometry and moving boundaries at different operating conditions can be confirmed. The SAS model is capable of resolving a major part of the incylinder turbulence intensity and is able to dynamically adjust the amount of resolved turbulence based on the local flow field and the spatial and temporal resolution. This benefits this model for industrial applications in comparison with classical LES models due to the possibility of using coarser meshes and moderate time step sizes. The overall agreement with experimental data measured at an optical accessible transparent engine regarding cycle averaged flow fields and rms fluctuations is very good. While there are some discrepancies in the intake phase, the matching of simulation and experiment at compression phase is very good. The motion of the characteristic tumble vortex was captured using the [[GAMMA].sub.1]-Criterion. The overall motion can be predicted by the numerical model.
Coupling a spray model, which was presented and validated with measurements at an optically accessible pressure chamber, with the cold flow model shows, that the mixture formation can be simulated. An analyzation of probability density function of the equivalence ratio showed, that there are remarkable differences in homogeneity of the mixture at right before and directly at TDC. Differences occur primarily on the lean side. High amounts of lean cells cannot be reduced throughout the tumble vortex breakdown.
Fluctuations of the equivalence ratio, temperature and NQVF at the spark plug at TDC can be resolved. Different location were analyzed and used to quantify the impact of fluctuations in the mixing process on the early stage of ignition using calculated values of the laminar flame speed of iso-octane. The COV of the laminar flame speed is nearly identical to the COV of the local equivalence ratio. The inhomogeneous temperature distribution can be neglected, since the impact on laminar flame speed is very small. Resulting variations of the laminar flame speed in the vicinity of the spark plug are expected to exceed 30% for the presented operating point. This will obviously have a direct impact on combustion.
The presented methodology shows, that SAS in combination with reasonable analyzation methods can be a useful and handy alternative to LES to analyze the incylinder mixing process with respect to cycle-to-cycle variations. This is of major importance for the industrial application of scale resolving CFD simulation. If sufficient computational resources are available, combining this model with optimization methods might be helpful to minimize CCVs. To go further and describe the whole process, accurate ignition and combustion models are needed in the computation.
[1.] Theile, M., Hassel, E., Thevenin, D., Buchholz, B., "LES of the cold flow of a DISI-Engine and validation with high-speed PIV measurements", Presentation at LES for Internal Combustion Engine Flows [LES4ICE], December 2014
[2.] Heywood, J.B., "Internal Combustion Engine Fundamentals", McGraw and Hill Series in Mechanical Engineering, McGraw-Hill, New-York, 1988
[3.] Dec, J. E., "Advanced compression-ignition engines - understanding the incylinder processes", Proceedings of the Combustion Institute 32:2727-2742, 2009
[4.] Granet, V., et. al., "Large-Eddy Simulation and experimental study of cycle-to-cycle variations of stable and unstable operating points in a spark ignition engine", Combustion and Flame 159(2012): 1562-1575
[5.] Johansson, B., "Cycle to Cycle Variations in S.I. Engines - The Effects of Fluid Flow and Gas Composition in the Vicinity of the Spark Plug on Early Combustion," SAE Technical Paper 962084, 1996, doi:10.4271/962084.
[6.] Brussies, E., "Simulation der Zylinderinnenstromung eines Zweiventil-Dieselmotors mit einem skalenauflosenden Turbulenzmodell", Ph.D Thesis, Mechanical Engineering Departement, Technical University Darmstadt, 2013
[7.] Enaux, B. et al. "LES study of cycle-to-cycle variations in a spark ignition engine", Proceedings of the Combustion Institute 33(2):3115-3122, 2011, doi:10.1016/j.proci.2010.07.038
[8.] Goryntsev, D. "Large Eddy Simulation of the Flow and Mixing Field in an Internal Combustion Engine", Ph.D. thesis, Mechanical Engineering Departement, Technical University Darmstadt, 2007
[9.] Haworth, D. C., Jansen, K. "Large Eddy Simulation on unstructured deforming meshes: towards reciprocating IC engines", Computer & Fluids 29(5):493-524, 2000, doi:10.1016/S0045-7930(99)00015-8
[10.] Rutland, C. J.: "Large-eddy simulations for internal combustion engines - a review:, International Journal of Engine Research 12(5):421-451, 2011, DOI: 10.1177/1468087411407248
[11.] Hasse, C., Sohm, V., Durst, B., ,,Numerical investigation of cyclic variations in gasoline engines using a hybrid URANS/LES modeling approach", Computers & Fluids 39:25-48, 2010
[12.] Imberdis, O., Hartmann, M., Bensler, H., Kapitza, L. et al., "A Numerical and Experimental Investigation of a DISI-Engine Intake Port Generated Turbulent Flow," SAE Technical Paper 2007-01-4047, 2007, doi:10.4271/2007-01-4047.
[13.] Bottone, F., Kronenburg, A., Gosman, D., Marquis, A., ,,Large Eddy Simulation of Diesel In-cylinder Flow", Flow Turbulence and Combustion 88:233-253, 2012, DOI: 10.1007/s10494-011-9376-6
[14.] Chang, C.-Y., "Development and Validation of Scale-resolving Computational Models Relevant to IC-engine Flow Configurations", Ph.D. thesis, Mechanical Engineering Department, Technical University Darmstadt, 2015
[15.] Basara, B., Poredos, A., and Gorensek, P., "Scale-Resolving Simulations of the Flow in Intake Port Geometries," SAE Technical Paper 2016-01-0589, 2016, doi:10.4271/2016-01-0589.
[16.] Basara, B., Krajnovic, S., Girimaji, S., Pavlovic, Z., "Near-Wall Formulation of the Partially Averaged Navier-Stokes Turbulence Model", American Institute of Aeronautics and Astronautics 49(12):2627-2626, 2011
[17.] Piscaglia, F., Montorfano, A., and Onorati, A., "A Scale Adaptive Filtering Technique for Turbulence Modeling of Unsteady Flows in IC Engines," SAE Int. J. Engines 8(2):426-436, 2015, doi:10.4271/2015-01-0395.
[18.] Menter, F.R., Kuntz, M., Bender, R. ,,A Scale-Adaptive Simulation Model for Turbulent Flow Predictions", American Institute of Aeronautics and Astronautics, 2003-767, DOI: 10.2514/6.2003-767
[19.] Menter, F.R., Egorov, Y. "The Scale-Adaptive Simulation Method for Unsteady Turbulent Flow Prediction. Part 1: Theory and Model Description", Flow, Turbulence and Combustion 85(1):113-138, 2010, DOI: 10.1007/sl0494-010-9264-5
[20.] Egorov, Y., Menter, F.R., Lechner, R., Cokljat, D. "The Scale-Adaptive Simulation Method for Unsteady Turbulent Flow Prediction. Part 2: Application to Complex Flows", Flow, Turbulence and Combustion 85:139-165, 2010, DOI: 10.1007/sl0494-010-9265-4
[21.] Buhl, S., Hartmann, F., Hasse, C.: Identification of Large-Scale Structure Fluctuations in IC Engines using POD-Based Conditional Averaging, Oil & Gas Science and Technology - Rev. IFP Energie Nouvelles, DOI: 10.2516/ogst/2015022
[22.] Menter, F.R., ,,Zonal Two Equation k-[omega] Turbulence Models for Aerodynamic Flows", American Institute of Aeronautics and Astronautics, 02/1993, DOI: 10.2514/6.1993-2906
[23.] Sommerfeld, M., "Computational Fluid Dynamics of Dispersed Multi-Phase Flows", (ERCOFTAC, 2008), ISBN 978-91-633-3564-8
[24.] Baumgarten, C., "Mixture Formation in Internal Combustion Engines, (Springer-Verlag Berlin Heidelberg New York, 2006), ISBN 978-3-540-30835-5
[25.] Beale, J. C., Reitz, R. D., "Modeling Spray Atomization with the Kelvin-Helmholtz/Rayleigh-Taylor Hybrid Model". Atomization and Sprays, 9:623-650, 1999, DOI: 10.1615/AtomizSpr.v9.i6.40
[26.] Fath, G.M., "Spray Atomization and Combustion", American Institute of Aeronautics and Astronautics, 1986-0136
[27.] Bader, A., Keller, P., Hasse, C., ,,The influence of non-ideal vapor-liquid equilibrium on the evaporation of ethanol/iso-octane droplets", International Journal of Heat and Mass Transfer 64:547-558, 2013, DOI: 10.1016/j.ijheatmasstransfer.2013.04.056
[28.] Lide David R., ed., CRC Handbook of Chemistry and Physics, Internet Version 2005, http://www.hbcpnetbase.com, CRC Press, Boca Raton, FL, accessed Jan. 2015
[29.] Kumzerova, E., Esch, T., Menter, F.: "Spray Simulations: Application of Various Droplet Breakup Models", presented at 6th International Conference on Multiphase Flow, Germany, July 9-13, 2007
[30.] Foucher, F., Landry, L., Halter, F., Mounaim-Rousselle, C., "Turbulent flow fields analysis of a Spark-Ignition engine as function of the boosted pressure", presented at 14th International Symposium of Laser Techniques to Fluid Mechanics, Portugal, July 7-10, 2008
[31.] Graftieaux, L.: "Combining PIV, POD and vortex identification algorithms for the study of unsteady turbulent swirling flows", Measurement Science and Technology 12(9):1422-1429, 2001
[32.] Janas, P., Wlkoas, I., Bohm, B., Kempf, A., ,,On the Evolution oft he Flow Field in a Spark Ignition Engine", Flow, Turbulence and Cobustion 1-18, 2016, DOI: 10.007/s10494-016-9744-3
[33.] Aleiferis, P.G., v. Romunde, Z.R.: "An analysis of spray development with iso-octane, n-pentane, gasoline, ethanol and n-butanol from a multi-hole injector under hot fuel conditions", Fuel 105:143-168, 2013
[34.] Cundy, M., Schucht, T., Thiele, O., Sick, V. ,,High-Speed laser-induced fluorescence and spark plug absorption sensor diagnostics for mixing and combustion studies in engines", Applied Optics, 48(4):B94-B104, 2009
[35.] Curran, H. J., Gaffuri, P., Pitz, W. J., Westbrook, C. K., "A Comprehensive Modeling Study of iso-Octane Oxidation", Combustion and Flame 129(3):253-280, (2002), DOI: 10.1016/S0010-2180(01)00373-X
FVTR GmbH, Germany
Ph. +49 381 498 9501
Fax: +49 498 9402
BDC - Bottom dead center
CAD - Crank angle degree
CARB - California Air Resources Board
CCV - Cyclic combustion variability
CFD - Computational Fluid Dynamics
DISI - Direct injection spark ignition
DPM - Discrete Phase Model
ER - Equivalence ratio
EVO - Exhaust valve opening
IVO - Intake valve opening
KHRT - Kelvin-Helmholtz Rayleigh-Taylor
LES - Large-Eddy Simulation
LEV - Low Emission Vehicle
LvK - Von Karman Length
OP - Operating point
PDA - Phase Doppler Anemometry
PDF - Probability density function
PIV - Particle image velocimetry
POD - Proper Orthogonal Decomposition
RANS - Reynolds-averaged Navier-Stokes
SST - Shear Stress Transport
TDC - Top dead center
TKE - Turbulent kinetic energy
[[GAMMA].sub.1] - Gamma 1 criterion
Martin Theile FVTR GmbH
Egon Hassel University of Rostock
Dominique Thevenin University of Magdeburg
Bert Buchholz University of Rostock
Karsten Michels and Martin Hofer Volkswagen AG
Table 1. Datasheet of the research engine used for CFD simulations Name value Stroke 80 mm Bore 74.5 mm Connecting rod 144 mm Compression ratio 7:1 Number of valves 4 Exhaust valve open 18.4[degrees] BBDC@ 0.2 mm lift Exhaust valve close 43.2[degrees] ATDC @ 0.2 mm lift Inlet valve open 42.1[degrees] BTDC @ 0.2 mm lift Inlet valve close 13.5[degrees] ABDC @ 0.2 mm lift Table 2. Description of the simulated operating points Operating Engine speed [P.sub.intake] [COV.sub.intake /exhaust] Point / 1/min /bar /% OP 1 1000 2 0.13/0.21 OP 2 2000 1 0.09/0.11 OP 3 3000 2 0.14/0.48 Operating No. Point cycles OP 1 10 OP 2 30 OP 3 10 Table 3. Operating points for spray injection in a constant volume chamber Operating [P.sub.chamber] [T.sub.chamber] [P.sub.fuel] Point /bar /K /bar [OP.sub.spray] 1 2.10 373 200 [OP.sub.spray] 2 0.62 373 200 [OP.sub.spray] 3 0.85 373 200 Operating [T.sub.fuel] [P.sub.gas] Point /K /kg/[m.sup.3] [OP.sub.spray] 1 363 1.961 [OP.sub.spray] 2 363 0.579 [OP.sub.spray] 3 363 0.794 Table 4. Description of the simulated operating points for analysis of mixture formation No. mesh type cell size time step size / [mm] /[[micro]s] 1 hexahedral 1 5 2 hexahedral 0.6 5 3 hexahedral 0.6 2 4 tetrahedral 0.6 2 Table 5. Description of the simulated operating points for analysis of mixture formation Name Engine speed [P.sub.intake] [m.sub.fuel]/ /1/min /bar mg [OP.sub.mixture] 1 1000 2 60.4 [OP.sub.mixture] 2 2000 1 34.2 [OP.sub.mixture] 3 3000 2 63.9 Name EOI/CAD No. b TDC cycle [OP.sub.mixture] 1 239 10 [OP.sub.mixture] 2 229 30 [OP.sub.mixture] 3 156 10 Table 6. Coefficient of variance of the equivalence ratio, temperature and NQVF shown in figure 24 Sphere radius [COV.sub.er] [COV.sub.T] [COV.sub.nqvf] r = 1 mm 7.91 0.73 12.43 r = 10 mm 2.66 0.23 12.47 Sphere radius [COV.sub.lamFS] r = 1 mm 7.96 r = 10 mm 2.70
|Printer friendly Cite/link Email Feedback|
|Author:||Theile, Martin; Hassel, Egon; Thevenin, Dominique; Buchholz, Bert; Michels, Karsten; Hofer, Martin|
|Publication:||SAE International Journal of Engines|
|Date:||Dec 1, 2016|
|Previous Article:||Evaluation of electrostatic screen battery for emissions control (ESBEC) with diesel emissions.|
|Next Article:||Studies on the effect of in-cylinder charge stratifications on high load HCCI combustion.|