© Y. Shekhawat et al., published by IFP Energies nouvelles, 2017

Introduction

Cycle-to-Cycle Variations (CCV) of flow and combustion in spark-ignition engines have been the subject of extensive research over the last few decades [1]. High-speed/high-resolution optical diagnostics and Large-Eddy Simulation (LES) have been brought to bear to understand the root causes of CCV, and significant progress has been made in that regard. Multiple sources of variability have been identified. These include local variations in the flow, temperature and mixture composition, and in the spark discharge, as well as global variations in equivalence ratio, dilution and trapped mass [2–5]. While these contributing factors have been recognized for some time, their interactions and influences on combustion remained unclear due to a lack of multi-diagnostic data and the difficulty of performing accurate LES simulations over a sufficiently large number of cycles. Lacour and Pera [5] and Baum et al. [6] performed multi-diagnostic experiments that allow a deeper understanding of the coupled physics in engines. And recently these interactions have been investigated using multi-parameter experimental [6,7] and simulation [6,8] approaches.

Earlier experimental and simulation studies (e.g., [9,10]) have shown that what happens during the earliest stages of combustion can determine, to a large extent, the subsequent combustion process for that cycle. Under realistic engine operating conditions, the velocity magnitude [7,8] and velocity gradient parameters [7] in the vicinity of the spark plug at the time of ignition have been found to correlate with the Indicated Mean Effective Pressure (IMEP) for the cycle. Motivated by the difficulties of making accurate pressure-based combustion measurements during the early flame development period, some LES studies have focused on early-flame-kernel development, toward developing an understanding of the governing factors that result in cyclic variation of flame growth. Granet et al. [11] demonstrated that initial flame convection within the spark plug gap is a probable reason for incomplete combustion cycles due to local quenching, while Goryntsev et al. [12,13] emphasized the superposition of flow variation and mixture quality for direct-injection engines, whose combined effect on combustion CCV was elucidated by LES.

For these reasons, the focus here is on the early burn: ignition through fully developed turbulent flame. A specific goal is to assess the predictive capability of two different ignition and turbulent flame propagation models, using the same CFD code and computational mesh. The approach is to compare the simulated and measured combustion using a single set of metrics that are applicable to both the simulated and measured data, focusing on the early burn. The metrics are derived from in-cylinder pressure and high-speed optical measurements, and are applied to perform equivalent analysis of the simulation data. The study was conducted in a research engine that exhibits high cycle-to-cycle variability. As will be shown later, fast and slow burning cycles are established during the Early Flame Development (EFD), defined here as the period between the Start Of Ignition (SOIgn) and a fully developed turbulent flame. Variation of the EFD duration (either between the experiments and simulations, or between individual cycles in either the experiments or the simulations) changes the combustion phasing, and thus is effectively the same as variability in the ignition timing. Thus, the turbulent flame propagation after the EFD encounters earlier or later thermodynamic states and turbulence properties from one cycle to another. In turn, the resulting combustion-phase variability affects the production of work (IMEP) and the time for completion of combustion (combustion efficiency), as revealed by experimental pressure-based combustion analysis [14,15]. Thus, improper simulation of EFD will result in improper pressure or mass burned fraction evolution, even if the turbulent flame propagation model is perfect. This motivates the focus on EFD for this study.

Historically, the EFD period has been characterized using the Mass-Burned Fraction (MBF) of the combustion progress, computed from an Apparent Heat Release (AHR) analysis based on the measured in-cylinder pressure trace. Typically, a crank angle corresponding to between one and ten percent MBF (CA01–CA10) is the earliest value that is sufficiently resolved experimentally, and that has been used as the metric to characterize the EFD variability. Although values as low as CA01 can be calculated from the pressure data, single-cycle analysis has high uncertainty due to recording noise and the large dynamic range requirements for the pressure transducer and recording system (here 40–2000 kPa). As a second EFD metric, the transition between slow and rapid flame growth has been measured optically. Beretta et al. [16] observed the early-flame-growth period between the start of ignition and the transition to a turbulent flame, which they defined in terms of the transition from a laminar to a turbulent characteristic time scale. This concept was introduced for modeling ignition in RANS simulations by Abraham et al. [17]. The duration of this laminar-to-turbulent transition period, τ_lam-turb, is measured here from the two-dimensional area projection of the three-dimensional flame volume, as computed from both the optical flame measurements and the LES. The laminar-to-turbulent transition is quantified as the time-scale transition revealed in the log(A_Flame) versus log(t) plots (where A_Flame is the time-varying flame area), as suggested by Arpaci et al. [18]. It is demonstrated that τ_lam-turb occurs at a MBF much less than one percent in the experiments.

An emphasis in this paper is to apply equivalent metrics to compare results from the experimental measurements with results from two different LES ignition and turbulent combustion models, focusing on the EFD period. As will be shown later, both models were found to produce higher variability than the experiment. Both the measured and simulated cycles show that the combustion phasing of early- and late-burning cycles is established by CA10. The projected flame-area growth rate between the start of the spark discharge (SOIgn) and CA10 is compared to τ_lam-turb, to CA01, and to the spark-plasma duration to estimate the time and length scales of the flame kernel during these periods. To assess the causes of the reasons for the disparity between the experiments and simulations, both the thermodynamic and the flow properties at SOIgn are compared.

Contributions of this paper are: (1) the first comparisons between combustion LES and experiment for a configuration that is fully accessible to other research groups; (2) comparisons of results from two different LES turbulent combustion models using the same mesh and numerics; and (3) combining pressure-based combustion analysis, optical diagnostics and LES to generate new physical insight into the early stages of combustion and the limitations of current modeling approaches. Specific questions that are addressed include:

• Are there specific local or global conditions that result in a fast cycle versus a slow cycle?

• Can useful criteria/metrics be established for the transition between ignition and fully developed turbulent flame propagation?

• What controls the heat release rate during the early burn?

• How well do current LES models represent this phase of combustion?

• What experiments/models are needed to make further progress in understanding and modeling the early burn?

1 Engine Experiments

The University of Michigan transparent combustion chamber (TCC-III) optical engine is illustrated in Figure 1 [19–21]. This is a single-cylinder, spark-ignited, two-valve, four-stroke-cycle, pancake-shaped-combustion-chamber engine with a transparent liner and piston. The geometric compression ratio is 10:1, and bore × stroke are 92 × 86 mm. This is a research engine that has been developed specifically to facilitate optical access for Particle-Image Velocimetry (PIV) and other high-fidelity laser-based diagnostics, as well as simple mesh generation for CFD. This geometrically simple engine configuration is known to exhibit highly variable in-cylinder flow. The combination of a pancake-shaped combustion chamber, large piston-to-valve diameter ratio (92 mm/30 mm), and undirected ports produces large CCV of the large-scale, motored flow structures [22,23]. Details of the geometric configuration and operating conditions, and extensive motored engine PIV data, are available at [21]. Several groups have been using the motored TCC engine data to test and validate physical models and numerical algorithms for in-cylinder LES [19,20,23–25].

The measured combustion data are from the fired TCC-III engine experiments conducted by Schiffmann [7]. There, the measured data for a range of operating conditions (lean, rich, N₂ dilution, methane fuel, propane fuel) were analyzed to explore the sources of CCV. This first comparison between simulated and measured combustion data is conducted at a robust combustion condition: 40 kPa Manifold Absolute Pressure (MAP), 1300 rpm, stoichiometric premixed propane-air reactants, and no N₂ dilution. Ignition was at 342 crank angle degrees after Top-Dead-Center exhaust (°aTDCe). The 40 kPa/101.5 kPa intake/exhaust manifold absolute pressures result in back-flow during exhaust-valve opening and valve overlap periods. Details of the engine and optical methods can be found in [7]. Engine and operating parameters key to this paper are provided in Table 1, and details are available at [21].

In-cylinder pressure measurements were recorded every 0.5 Crank-Angle Degrees (CAD) of rotation. The AHR analysis was performed on both the measured and simulated pressure traces using the same crank-angle resolution and similar values of the polytropic coefficients n; for the experiments, individual-cycle values of n were used [7], whereas averaged values were used for the LES. The instantaneous heat release for crank-angle interval i is given by (1)

The cumulative heat release (Q) is then computed by integrating equation (1), and the AHR-based mass-burned fraction is given by MBF = Q/Q_total, where Q_total = Q_max.

Imaging of OH* chemiluminescence (Fig. 1) was used to measure the area of the two-dimensional projection of early flame growth with 0.5 × 0.5 mm spatial resolution. High-speed imaging was used to acquire images every 2 CAD from SOIgn (342 °aTDCe) to top-dead-center combustion (corresponding to an average MBF of approximately CA10), at which time the flames had grown outside of the field of view. The OH* signal was filtered by a combination of optical filters (UG11 and a 330 nm short-pass filter) that is centered on 307 nm with a 25 nm wide 10% cut-off to isolate the 306.4 nm band from emissions of other major combustion products over 330 nm [26]. To account for spatial differences in intensifier efficiencies, all acquired images were first white-field corrected and a spatial calibration was applied. Then a binarization threshold was determined by the average value of the maximum values in each pixel column. The parameterization of this threshold improves capturing the burned-gas-area region and allows a higher degree of automation. The flame-kernel growth rate can be quantified in this manner prior to the time when acceptable pressure measurements are available, and the measured growth rates then can be compared with the simulated flame-kernel growth rates.

High-speed PIV was used to measure the two-dimensional, two-component velocity in the plane y = −4.5 mm near the spark plug (Fig. 1). The present analysis focused on data from a 6 × 12 mm region (x:z = −2 → −8 mm: −8 → +4 mm), as shown with the white box in Figure 1. The spatial resolution was 1 mm, and was determined by the 1 mm thick laser sheet and the 1 × 1 mm (32 × 32 pixel) interrogation spots on a 0.5 × 0.5 mm grid (50% interrogation-spot overlap). The velocity dynamic range was 0.6–25 m/s, as determined by the 0.2-to-8-pixel PIV resolution. The velocity gradients were computed using central differencing on the grid interior, with one-sided differencing at the boundaries. The measurement plane was offset by 4.5 mm from the spark gap (y = −4.5 mm, Fig. 1) to avoid scattering from the spark-plug elements. Details of the PIV measurement procedures are provided in [7].

Figure 1 Schematic illustration of the TCC-III engine showing the coordinate system, PIV measurement plane (y = −4.5 mm) and the flame imaging fields-of-view. The white box in the PIV field-of-view image shows a 12 × 6 mm region near the spark plug (x = 4 to −8 mm, z = −2 to −8 mm) in which flow parameters have been averaged in the subsequent analysis.

2 Numerical Simulations

LES has been performed using STAR-CD v4 [27]. The computational domain includes the in-cylinder region, intake and exhaust ports, and intake and exhaust plenums. The average in-cylinder mesh size is approximately 0.75 mm, with finer meshes near the spark-plug tip and in the valve-curtain regions. Time-dependent pressure and temperature boundary conditions from a one-dimensional model are prescribed at the plenum inlet/outlet, with no CCV [19]. A standard Smagorinsky subfilter-scale turbulence model is used. Results from two different ignition and turbulent flame propagation models (designated as Model 1 and Model 2) are presented, with all other numerical, flow and heat transfer modeling aspects being the same. For each combustion model, 60 engine cycles have been simulated. For Model 1, three blocks of 20 consecutive cycles were simulated, whereas for Model 2, four blocks of 15 consecutive cycles were simulated. The same initial conditions were used for both models, which were obtained from different cycles of a preliminary LES study where a different subfilter-scale turbulence model and turbulent combustion model were used. No patterns of correlation between different blocks of cycles were evident, as shown in Figure 2.

Figure 2 Maximum in-cylinder pressure versus cycle number for Model 1 (top: 3 sets of 20 consecutive cycles) and Model 2 (bottom: 4 sets of 15 consecutive cycles).

2.1 Model 1

Model 1 uses a simple energy deposition model for ignition, and a modified Thickened Flame Model (TFM) [28] implemented by the authors for turbulent flame propagation. The energy deposition model creates a spherical temperature source with a radius of 0.75 mm, a temperature of 3000 K, and a duration of 4.0 crank angle degree. This high-temperature source initiates a flame based on a two-step global propane-air chemical mechanism [11], which has been tuned to give the correct laminar flame speed and adiabatic flame temperature for the engine conditions of interest. After the flame initiation, the thickened flame model takes over to propagate the flame. The essential features of the thickened flame model can be understood in terms of classical thermal theory for laminar premixed flames, according to which the laminar flame speed is proportional to the square root of the product of the molecular diffusivity D and the reaction rate ω, and the laminar flame thickness is proportional to the square root of the quotient of D and ω: (2)

The first essential idea in TFM is to multiply the diffusivity D and to divide the chemical source term ω in the species and enthalpy equations by the same factor F. From equation (2), it follows that this increases the flame thickness by the factor F (for F > 1), while the laminar flame speed (or equivalently, the mass rate of reactant consumption) remains unchanged. This suggests an approach by which one can recover the correct burn rate in an “under-resolved” numerical simulation of a flame where the underlying flame structure is not fully resolved: choose F to be the ratio of the actual mesh size to the mesh size that is needed to fully resolve the underlying flame structure (typically ∼0.1 ). In LES of turbulent flames, an additional factor is needed to account for the wrinkling of the flame (thereby enhancing the local burn rate) by unresolved turbulence scales. The second essential element in TFM, then, is to maintain the correct overall burn rate (or turbulent flame speed) in a numerical simulation of a turbulent flame, even when all turbulence scales are not fully resolved. This effect is implemented through a multiplicative factor Ξ that is applied to both the diffusive term and the chemical source term in the species and enthalpy equations. Multiplying both D and ω in equation (2) by the same factor Ξ increases the flame speed (or mass rate of reactant consumption) by the factor Ξ (for Ξ > 1), without changing the flame thickness. In contrast to F, Ξ has a rather complicated functional form that depends on local turbulence parameters, and that is informed by theoretical analysis and by direct numerical simulations. The resulting species mass fraction (Y_m) equation is then: (3) and similarly for enthalpy. The factors F and Ξ are applied only within the turbulent flame brush. Several models for Ξ have been formulated. Charlette et al. [29] proposed the following power-law: (4)

Here is the combustion filter size, β is a model constant, the function Γ represents the integration of turbulent strain rate induced by scales smaller than Δ (which is fitted by a function [29]), and is the subfilter-scale velocity fluctuation. An expression for the latter was proposed by Colin and Truffin [30]: (5) where  is a model constant. The resolved velocity gradients are felt by the turbulent combustion model through the Smagorinsky subfilter-scale turbulence model and through . Two important modifications for applications of TFM in engines are to limit the turbulence length scale in vicinity of walls (the length scale cannot be larger than the distance to nearest wall), and to disable the unrealistically high rate of wall heat transfer that is otherwise implied by the model when the flame reaches a wall. Details can be found in [31].

2.2 Model 2

Model 2 uses the Imposed Stretch Spark-Ignition Model (ISSIM-LES) [30], implemented in STAR-CD by the authors and coupled with the ECFM-LES turbulent combustion model [32]. The model is based on the theoretical ability of the flame surface density (FSD) approach to account for flame kernel growth to a fully developed flame, from the first appearance of a minimum self-sustained flame nucleus. The simplified model of the coil ignition system from Duclos and Colin [33] is adopted, accounting for electric energy discharge after the breakdown electrode voltage has been exceeded. The electric circuit of the TCC-III engine delivers 36.5 mJ per cycle, and this value is used in the calculations. The model reproduces the spark-to-gas energy transfer rate in the glow phase as a function of voltage and electric intensity, so that a finite-duration (non-instantaneous) ignition process is simulated. Energy losses to the electrodes lead to an energy transfer efficiency of approximately 60%, following Verhoeven's measurements [34]. Whereas the model for the electric circuit is able to account for re-striking events, no information on this is available from the experiments. The spark plug is fully meshed in the 3D model, and the energy at breakdown is deposited at the midpoint of the electrode gap (1.15 mm). A corresponding amount of burned gas is initialized, equal to the mass of a cylinder with radius 2δ_l, the electrode gap as the initial length and unburnt density ρ_u. The modeling of a reduced size flame kernel is a challenging task for the FSD-LES equation. For this purpose, a modified FSD-LES equation was proposed in [32]. With r_b denoting the burned-gas (flame) radius and  the flame surface density, the term  accounts for subfilter-scale stretch during ignition. As stated in [30], the hypothesis that the wrinkling factor Ξ = 1 is made, given the low turbulence intensity at ignition time. The modified FSD-LES equation is: (6)

The transition from subfilter-scale-only to subfilter-scale and resolved-scale contributions to flame area corresponds to the continuous transition from ignition to a fully propagating turbulent flame. The use of a unique transport equation for FSD from the early stages allows the simulation of ignition events occurring in regions of high flow velocity and/or velocity gradients. The transition is governed in the model by an algebraic α function, which depends on the filtered flame radius and is expressed as in equation (7). (7) The flame radius is calculated from the equivalent enflamed spherical volume V_b. This is defined as shown in equation (8); it is function of the mass-based progress variable , where  and  are the Favre-averaged concentrations of fuel and fuel tracer: (8)

At ignition (a ≅ 0), the ignition source term is active and accounts for flame surface production due to subfilter-scale stretch. As the enflamed volume grows and α increases from 0 to 1, the ignition source term is gradually suppressed, while the C_sgs, , C_res and S_res terms are activated. These are the subfilter-scale and resolved-scale curvature, the flame surface density production and the resolved strain contributions to the FSD equation. These terms restore the standard flame-surface-density-based LES turbulent flame propagation equation. The use of a transport equation for the entire flame history starting from ignition allows for full interaction between the local flow conditions (resolved velocity and velocity gradients) and the initial flame growth. The gradient of the resolved velocity affects the resolved strain , while the subfilter-scale velocity fluctuation u′ is accounted for through the efficiency function  that was proposed in [28], and which is included in the subfilter-scale strain term S_sgs. In equation (9), Δ is the filter length, and s_l and δ_l are the laminar flame speed and thickness, respectively: (9)

Finally, a region-specific wall temperature of 1000 K is imposed for the spark plug metal electrodes, given the impact of heat transfer from the initially small enflamed kernel to the surrounding solid parts.

3 Results

The comparisons between experiments and simulations are organized into four subsections: analysis of thermodynamic conditions at time of ignition, analysis based on apparent heat release and related quantities, analysis based on experimental flame imaging and analogous quantities extracted from the simulations, and analysis of the velocity fields in the vicinity of the spark-plug gap at the time of ignition.

3.1 Thermodynamic Analysis

The first comparison between experiments and simulations is the thermodynamic state of the mixture at SOIgn, 342 °aTDCe. The log P−log V plots of Figure 3 show 60 cycles of experimental data and 60 cycles from each model through compression, combustion and expansion. This provides a comparison of the volume-averaged thermodynamic states experienced in the experiments versus those experienced by the two models. Looking first at the compression stroke, it can be seen that there is imperceptible CCV in both the measurements and the simulations, and that the measured data suggest a larger value for the polytropic coefficient n compared to the models: n_Meas = 1.33, n_Mod1 = 1.31, n_Mod2 = 1.31. There is very little droop in the expected straight lines, as indicated by the very high values of the measured and simulated R². The lower polytropic coefficients for the LES can only be caused by different trapped mass or by heat loss, since there is no blow-by in the LES and the air/fuel composition is well known. The unknown temperature of the quartz cylinder in the experiments and uncertain wall-heat-transfer model in the LES are the most likely reasons for the disagreement between the measured and simulated polytropic coefficient values. It can be estimated from the polytropic compression that this difference in the values of n would result in a 6% lower bulk temperature for the LES compared to experiment (out of approximately 600 K) at SOIgn = 342 °aTDCe.

The pressure CCV near TDC are also revealed in Figure 3. From these, it can be seen that some of the LES cycles burn quite late, and that in turn yields higher CCV in the expansion stroke in the simulations compared to the experiments.

Figure 3 log P versus log V diagrams of the compression, combustion and expansion strokes for the measured data (top), Model 1 (center) and Model 2 (bottom).

3.2 Apparent Heat Release Analysis

The combustion variability is more evident in Figure 4, where measured and computed Apparent Heat Release Rate (AHRR) traces are compared. There, 60 measured cycles are shown, so that the sample size is the same as in the simulations. It will be shown later that these 60 measured cycles are representative of the 1500 total available experimental cycles (500 cycles each from three separate tests). The measured results in Figure 4 show that combustion has completed by 30 CAD after TDC for all cycles. This was suggested by Ayala and Heywood [4] to be a necessary condition for stable engine operation. Both LES models exhibit longer combustion durations and higher CCV compared to what is observed in the experiments. In particular, the simulations (especially Model 1) show some extremely slow cycles that are not seen experimentally for this operating condition.

It is customary to use the MBF obtained from the pressure-based cumulative AHR, to quantify the 10-, 50- and 90-percent MBF periods (CA10, CA50 and CA90). Figure 5a compares the individual-cycle MBF based on the cumulative AHR with the MBF based on the actual fuel consumption for LES Model 1, and Figure 5b shows the corresponding CA10, CA50 and CA90 values. Results for Model 2 are similar (not shown). Together, Figure 5a and b demonstrate that the pressure-based AHRR analysis of equation (1), which has typically been applied to experimentally measured data, is a precise metric that accurately reflects the actual fuel burn rate. It is therefore useful as a single metric that can be applied to both the measured and simulated data.

Figure 6 shows the range of CA10 CCV, and the correlations between CA50, CA90, and the Location of Peak Pressure (LPP) with CA10. Linear correlations of CA50, CA90 and LPP with CA10 are evident in both the experiments and the LES. The experiments exhibit a consistently smaller range of values for CA10 and CA90 (lower CCV) compared to the simulations. However, in all cases the burning after CA50 becomes more variable, as indicated by the lower correlation of CA10 with CA90 compared to CA50. Although the ranges of the values from the simulations overlap with the measured values, it is clear that both models exhibit some cycles with considerably slower burns compared to the measurements. More importantly for the present purposes, Figure 6 demonstrates that the CCV of combustion is highly correlated with the early flame kernel development time, even for the very-late-burning simulated cycles (evident from the R² values). Therefore, it is essential to correctly capture the early flame development period in the simulations, and to better understand what is happening during this critical period. The subsequent comparisons between experiments and simulations focus on the early flame development.

While it is possible to compute the AHRR from measured pressured data prior to CA10, experimental values of AHRR at very early times can be unreliable because of inaccuracies and low signal-to-noise ratios. Figure 7a shows this breakdown for the very early burn (CA005). Figure 7b shows that CA01 correlates linearly with CA10 in both the experiments and the simulations. As already noted, the experiments are difficult for very early times. And in the LES models used here, factors such as discharge repeatability and heat transfer to electrodes, are not fully taken into account; therefore, a correlation between CA01 and CA10 might be expected in the simulations, to some extent. The arbitrary imposition of a cycle-independent ignition energy might reduce the CA01/CA10 correlation, but possibly for the wrong reason. With these caveats, the results suggest that the phasing of CA10 (and from what has been shown earlier, the phasing of the entire combustion event) has largely been established by CA01 for this operating condition. Thus, it is of interest to explore additional measures to dissect the ignition and flame propagation processes prior to CA01. For this purpose, flame imaging was used to capture the early-flame growth, and to estimate the notional laminar-to-turbulent transition time, τ_lam-turb, suggested by Abraham et al. [17], as described in the following subsection.

Figure 4 Measured and computed individual-cycle AHRR versus crank angle. Top: measured. Center: Model 1. Bottom: Model 2.

Figure 5 Comparisons of MBFs based on the cumulative AHR versus on the fuel consumption rate, for Model 1.

Figure 6 Scatter plots of CA50, CA90 and Location of Peak Pressure (LPP) versus CA10. Left: Measured. Center: Model 1. Right: Model 2.

Figure 7 (a) Measured MBF at crank angles prior to CA10, demonstrating limitations of pressure-based AHRR analysis prior to CA01. (b) Comparisons of CA01 and CA10 distributions between simulated (Model 1 and Model 2) and measured (labeled “UM”) cycles.

3.3 Flame Imaging Analysis

The flame growth was quantified here as the projected 2D area of the 3D burned-gas volume, as measured from the OH* images and the simulations (Fig. 8a) every 2 CAD. For the projected flame area in the simulations, the 3D flame surface was defined as the location of 50% local mass burned fraction, and the 2D projected area was computed as in the experiments. This is a somewhat less-than-perfect comparison between experiments and simulations, because of the missing out-of-plane flame-growth information. Nevertheless, it is of interest to see whether or not this relatively simple approach reveals anything of interest. Applying the concept of Arpaci et al. [18], the transition time τ_lam-turb is defined here as the time between SOIgn and max(d²A_flm/dt²), as shown in example of Figure 8b; therefore the roll-off of the flame growth after 1 ms is due to the limited experimental field of view. Figure 8c shows that this transition occurs well before CA01.

It is recognized that the single view that is available experimentally is limited in its ability to capture out-of-plane growth, for example. However, this metric is intended only to identify the overall growth rate, not the absolute volume. This metric was found to correlate with CA10 and to discriminate between the stoichiometric and rich/lean/dilute mixtures of the experimental measurements of Schiffmann [7]. Figure 9 shows that there is a weak linear correlation between τ_lam-turb and CA10 (the bimodal distribution of growth rates is believed to be an artifact of the low sampling rate (2 CAD) compared to the growth rate). This may be due to multivariable influences as well as the single view of the measurement. Nonetheless, τ_lam-turb identifies the transition from one growth rate to another, and in principle can be applied equally to both the measured and simulated data. However, a sharp transition in growth rate was not observed in the simulations for either model (not shown).

To further dissect the EFD period, the measured cycle-to-cycle spark-duration distribution is compared to the measured duration of CA10 (CA00–CA10) in Figure 9. The plasma duration distribution is presented here as a reminder that both electrical and chemical energies contribute to the early flame kernel during this period. Figure 9 demonstrates that the glow-discharge ignition system used here creates a plasma whose duration is longer than τ_lam-turb, and that it lasts until nearly CA10. This suggests that the flame growth measured by the OH* emissions have commenced well before the plasma has ended. Furthermore, neither spark duration nor τ_lam-turb correlate strongly with CA10. This figure also illustrates that the 60 measured cycles that were selected to compare with the 60 computed cycles from each model are representative of the 1500 available measured cycles.

A length scale that is expected to be relevant for the early flame development is computed as the diameter of a circle that has the same projected 2D area as that of the flame kernel. Figure 10 shows the ensemble-averaged equivalent diameter of the flame as a function of time after start of ignition, for the experiments and the two LES models. The range of the CCV in this quantity (dotted lines) is quantified by the minimum and maximum areas observed in the three data sets. Also shown are the average values of τ_lam-turb, CA01, End of Spark (EoSpk), and CA10. From this figure, estimates of the relevant flame-kernel length- and time-scales during the period between SOIgn and CA10 can be deduced. Since τ_lam-turb has been shown to be poorly correlated with CA10, while CA01 is highly correlated with CA10, it can be deduced that the correlation with CA10 develops in the period between 0.4 ms and 1.5 ms (CA01), where the flame diameter is between 5 and 15 mm. This can be compared to the estimated integral time scales of 1–1.3 ms and integral length scales of L_ux, L_uz, L_wz and L_wx equal to 10–20 mm, 8 mm, 3 mm and 5 mm, respectively [7], where L_i,j denotes the correlation length for the ith velocity component in the jth direction. This then suggests that simulated and measured flow properties at these scales could be compared, to develop further insight into the early flame development.

Figure 8 (a) Examples of measured and simulated 2D projected flame kernel area. A single-cycle example of the burned-gas area evolution for a measured cycle is used to quantify the laminar-turbulent transition time, τlam-turb, shown by the red-diamond in (b), which occurs much earlier than the 1% MBF from the AHR analysis (c).

Figure 9 Measured spark-discharge duration and τlam-turb distributions from 1500 cycles (500 cycles from each of three independent tests − in black), and the 60 sampled cycles (in red) that are used to provide a sample size that is the same as that available from the simulations.

Figure 10 Ensemble-averaged (solid lines) and cycle-to-cycle range (dotted lines) of the flame's equivalent diameter for all measured (1500) and simulated (60 each) cycles. The vertical lines mark crank angles of interest, from which the temporal and spatial scales of the early flame-kernel growth can be deduced.

3.4 Local Flow at Spark Gap Analysis

Comparisons between computed and measured mean and rms velocity fields for motored and firing conditions are reported in [31], and the reader is referred there for detailed quantitative comparisons. An example of mean velocity fields on two cutting planes that pass close to the spark gap at an instant shortly before SOIgn are shown in Figure 11, to give a general impression of the large-scale flow structure. There are quantitative differences between the two LES models, between LES and PIV, and between the motored and firing pre-ignition flow fields, but the general flow patterns are similar.

Multiple earlier experimental and simulation studies have shown that local flow conditions in the vicinity of the spark gap at the time of ignition affect the early flame kernel development and subsequent combustion process [7–13]. In particular, flame plasmas and kernels can be advected away from the spark plug, which affects the combustion phasing. Thus, it is prudent to compare not only the thermodynamic state between the experiments and simulations at the time of ignition (Fig. 3), but also the velocity and velocity gradients. Scatter plots of velocity magnitude, |V|, and the direction of the velocity in the 2D measurement plane (Fig. 1), |u| and |w|, are shown in Figure 12. Both the measured and simulated values have been spatially averaged over the 12 mm × 6 mm region illustrated in Figure 1. Therefore, these represent low-pass spatially filtered values (over a spatial scale of ∼10 mm) in the vicinity of the spark plug. As noted earlier, the measurement plane does not pass directly through the spark-plug gap (it is offset by 4.5 mm). Nevertheless, the measured values of |V| were found to correlate with CA10 in the statistical analysis of Schiffmann [7]. This offset does not detract from the point of the comparison here, which is to determine whether or not the simulations are predicting similar velocities near the spark plug. Figure 12 shows that the range of |V| values from the simulations tend to lie within the lower half of the measured range of values. Furthermore, the direction of the large-scale flow (as indicated by the trend in the |w| versus |u| plot) is consistent between the experiment and simulations, although again, the overall lower velocity magnitudes in the simulations are evident.

Since Figure 12 demonstrates that there are differences between the measured and simulated large-scale velocity magnitudes near the spark plug, it is reasonable to assess if the measured and simulated resolved-scale velocity gradients are similar or significantly different. The resolved-scale velocity gradients play an important role in subfilter-scale turbulence and turbulent combustion modeling, as discussed earlier in the modeling section. In Figure 13, measured and computed distributions of the two-component resolved-scale velocity gradient components over the 12 mm × 6 mm measurement region displaced 4.5 mm from the spark-plug gap (the same region for which area-averaged velocities are shown in Fig. 12) are plotted. These distributions were sampled for all 60 simulated cycles and for the 60 representative measured cycles. The measured (PIV) velocity gradients were computed directly from 1 mm-resolution interrogations spots on a 0.5 mm grid, whereas the simulation distributions were sampled on the computational-grid and subsequently mapped onto the PIV grid. The most striking feature of Figure 13 is the significantly narrower distribution of measured velocity gradient values compared to the simulated velocity gradient values.

To assess if the larger gradients from the simulations are the result of higher computational-grid resolution compared to the experimental spatial resolution, a linear correlation was performed between the gradients mapped from the computational grid to the PIV grid and gradients computed from the remapped computed velocities. The linear correlation demonstrated that the gradients computed from the remapped velocity were only 8% lower than the higher-resolution simulation gradients, with a linear correlation coefficient of 0.98. Thus, it was concluded that the narrower distributions of the measured gradients are not due to the remapping. Rather, Figure 13 reveals the limited ability of PIV measurements to resolve steep gradients. This may be the result of three factors: the error between the particles' motion and the local gas-phase flow velocity; the limited dynamic range of the instantaneous velocity measurements (0.2 pixels out of 8 = 40×); and the difficulty of PIV-correlation-based peak-detection in regions of steep velocity gradients (vector drop out). Thus, based on these measurements, it is not possible to conclude if the simulations accurately predict the resolved-scale gradients in the flow. Nevertheless, since the two LES models yielded different combustion phasing (e.g., Fig. 4), it is of value to compare the computed gradient distributions between the two models.

To assess the velocity-gradient distributions that the LES models present to the ignition/combustion models, all nine velocity-gradient-tensor components were sampled and compared from the two simulations in a three-dimensional volume around the spark plug at the time of ignition. The distributions were sampled from all 60 cycles, on a 0.5 mm grid in a rectangular prism volume centered at the spark-plug gap: x = −8 mm → 4 mm, y = −4.5 mm → 4.5 mm, z = −8 mm → −2 mm. In addition, the distributions from the 10 fastest burning and 10 slowest burning cycles are plotted in Figure 14. These distributions are the fraction of the total sample in the histograms (bin size of 250 1/s) normalized by the number of samples. Thus, the fraction of the total sample (370k samples) can be compared to the 37k samples of the 10 fastest and slowest cycles, respectively. The distributions displayed in Figure 14 show that Model 2 has more low-magnitude gradients (−1000/s to +1000/s) compared to Model 1; Model 1 yields broader distributions with larger excursions from the mean. Figure 14 also shows slight variations in skewness between the two simulations. As demonstrated by Petersen and Ghandhi [35], the turbulent dissipation in engine flows is expected to occur in relatively localized high-shear regions between larger structures. Thus, the occurrence of the largest gradients from Eulerian measurements is expected to be infrequent compared to the smaller magnitude gradients at the interior of the larger structures. This is consistent with the distributions in Figure 14.

To better reveal the large dynamic range and relative frequency of occurrence of high-magnitude gradients, it is of value to replot the fractional frequency in Figure 14 using a logarithmic scale (Fig. 15). Figure 15 more clearly reveals the comparatively smaller range of velocity gradients presented to the ignition and turbulent combustion models by Model 2 compared to Model 1. The reason for this difference is unclear, at this point. It is perhaps more revealing to compare the fastest and slowest cycles from each model. Figures 14 and 15 show that for both models, the slowest cycles correspond to a narrower distribution of velocity gradient magnitudes (fewer samples with high velocity gradient magnitudes) compared to the mean distribution, while the fastest cycles correspond to a broader distribution of gradient magnitudes (more samples with high velocity gradient magnitudes) compared to the mean. Higher velocity gradients may lead to accelerated flame growth, and exposing the early kernel to a broader range of fluctuating velocity gradients may also contribute. Conversely, lower velocity gradient fluctuations and the absence of high velocity gradient magnitudes result in slower early flame growth.

Figure 11 Comparisons between measured and simulated mean velocity fields on two cutting planes at 340 °aTDCe (2° before SOIgn).

Figure 12 Comparisons between measured and simulated velocity magnitude |V| distributions with the duration of CA10 (left), and 2D velocity components |u| and |w| showing distribution of the velocity direction (right). In all cases, the velocities have been spatially averaged over a 12 mm × 6 mm region in the vicinity of the spark-plug gap (y = −4.5 mm, Fig. 1) at SOIgn (342 °aTDCe). For the experiments, all 1500 cycles and 60 representative cycles are shown, the latter to provide a clearer visual comparison with the 60-cycle LES results.

Figure 13 Measured and computed two-component velocity gradient distributions in a measurement region in the vicinity of the spark plug (Fig. 1) at SOIgn (342 °aTDCe). Distributions of all four 2D velocity gradient components are shown.

Figure 14 Distributions of nine components of the velocity gradient tensor in the volume centered at the spark plug at SOIgn (342 °aTDCe) for two different LES models. The distributions are sampled at all volume nodes in the region for the 60 simulated cycles (total sample size of 370k). Three distributions are shown for each component: all 60 cycles, the 10 fastest-burning cycles, and the 10 slowest-burning cycles.

Figure 15 Distributions of nine components of the velocity gradient tensor in the volume centered at the spark plug at SOIgn (342 °aTDCe) for two different LES models. The distributions are sampled at all volume nodes in the region for the 60 simulated cycles (total sample size of 370k). Three distributions are shown for each component: all 60 cycles, the 10 fastest-burning cycles, and the 10 slowest-burning cycles. This is the same as Figure 14, but plotted using a logarithmic scale.

Conclusion

Pressure-based combustion analysis, optical diagnostics and LES have been used to generate new physical insight into the early stages of combustion in homogeneous-charge spark-ignition engines. The emphasis has been on developing strategies for making quantitative comparisons between high-speed/high-resolution optical diagnostics and LES using common metrics for both the experiments and the simulations, focusing on the important early flame development period.

These are the first comparisons between combustion LES and experiments for the TCC engine, an engine configuration that is freely accessible to other research groups. Results from two different LES ignition turbulent flame propagation models have been presented, using the same numerical methods and computational mesh. The results reveal strengths and limitations of the experimental diagnostics and LES models. Specific findings are as follows:

• Local flow conditions in the vicinity of the spark gap at the time of ignition play a significant role in determining the subsequent combustion event for the cycle. This largely confirms results from multiple earlier experimental and simulation studies.

• Two widely used LES models yield rather different combustion behavior compared to experiments. Both models show somewhat lower velocity magnitudes in the vicinity of the spark gap at the time of ignition compared to experiment, and both models show higher combustion CCV compared to experiment. In particular, the models give some very slow-burning cycles that are not observed in the experiments. Nevertheless, the comparisons between models, and between models and experiments, provide useful insight into the early flame development and into the limitations of the models and the experiments.

• In the simulations, fast-burning cycles correspond to cycles where there is a broader distribution of velocity gradient magnitudes (more samples with high velocity gradient magnitudes) in the vicinity of the spark-plug gap at the time of ignition compared to the average distribution. Conversely, slow-burning cycles correspond to cycles where there is a narrower distribution of velocity gradient magnitudes in the vicinity of the spark-plug gap at the time of ignition compared to the average distribution. Higher velocity gradient magnitudes may lead to accelerated flame growth, and exposing the early kernel to a broader range of fluctuating gradients may also contribute. To explore this experimentally, it will be necessary to modify the PIV technique that has been used here to better resolve steep local velocity gradients.

• The fate of the combustion event for each cycle is largely determined by what happens between the laminar-to-turbulent transition time and CA01. This suggests a range of temporal and spatial scales over which subsequent comparisons between experiments and simulations should be performed.

Acknowledgments

Coauthors from the Pennsylvania State University and the University of Michigan gratefully acknowledge funding from GM Global R&D. Coauthors from the Pennsylvania State University and the University of Modena and Reggio Emilia gratefully acknowledge Siemens PLM Software for providing access to, and support for, their STAR-CD CFD code.

References

All Tables

All Figures

	Figure 1 Schematic illustration of the TCC-III engine showing the coordinate system, PIV measurement plane (y = −4.5 mm) and the flame imaging fields-of-view. The white box in the PIV field-of-view image shows a 12 × 6 mm region near the spark plug (x = 4 to −8 mm, z = −2 to −8 mm) in which flow parameters have been averaged in the subsequent analysis.
In the text

	Figure 2 Maximum in-cylinder pressure versus cycle number for Model 1 (top: 3 sets of 20 consecutive cycles) and Model 2 (bottom: 4 sets of 15 consecutive cycles).
In the text

	Figure 3 log P versus log V diagrams of the compression, combustion and expansion strokes for the measured data (top), Model 1 (center) and Model 2 (bottom).
In the text

	Figure 4 Measured and computed individual-cycle AHRR versus crank angle. Top: measured. Center: Model 1. Bottom: Model 2.
In the text

	Figure 5 Comparisons of MBFs based on the cumulative AHR versus on the fuel consumption rate, for Model 1.
In the text

	Figure 6 Scatter plots of CA50, CA90 and Location of Peak Pressure (LPP) versus CA10. Left: Measured. Center: Model 1. Right: Model 2.
In the text

	Figure 7 (a) Measured MBF at crank angles prior to CA10, demonstrating limitations of pressure-based AHRR analysis prior to CA01. (b) Comparisons of CA01 and CA10 distributions between simulated (Model 1 and Model 2) and measured (labeled “UM”) cycles.
In the text

	Figure 8 (a) Examples of measured and simulated 2D projected flame kernel area. A single-cycle example of the burned-gas area evolution for a measured cycle is used to quantify the laminar-turbulent transition time, τlam-turb, shown by the red-diamond in (b), which occurs much earlier than the 1% MBF from the AHR analysis (c).
In the text

	Figure 9 Measured spark-discharge duration and τlam-turb distributions from 1500 cycles (500 cycles from each of three independent tests − in black), and the 60 sampled cycles (in red) that are used to provide a sample size that is the same as that available from the simulations.
In the text

	Figure 10 Ensemble-averaged (solid lines) and cycle-to-cycle range (dotted lines) of the flame's equivalent diameter for all measured (1500) and simulated (60 each) cycles. The vertical lines mark crank angles of interest, from which the temporal and spatial scales of the early flame-kernel growth can be deduced.
In the text

	Figure 11 Comparisons between measured and simulated mean velocity fields on two cutting planes at 340 °aTDCe (2° before SOIgn).
In the text

Figure 12 Comparisons between measured and simulated velocity magnitude |V| distributions with the duration of CA10 (left), and 2D velocity components |u| and |w| showing distribution of the velocity direction (right). In all cases, the velocities have been spatially averaged over a 12 mm × 6 mm region in the vicinity of the spark-plug gap (y = −4.5 mm, Fig. 1) at SOIgn (342 °aTDCe). For the experiments, all 1500 cycles and 60 representative cycles are shown, the latter to provide a clearer visual comparison with the 60-cycle LES results.

In the text

	Figure 13 Measured and computed two-component velocity gradient distributions in a measurement region in the vicinity of the spark plug (Fig. 1) at SOIgn (342 °aTDCe). Distributions of all four 2D velocity gradient components are shown.
In the text

Figure 14 Distributions of nine components of the velocity gradient tensor in the volume centered at the spark plug at SOIgn (342 °aTDCe) for two different LES models. The distributions are sampled at all volume nodes in the region for the 60 simulated cycles (total sample size of 370k). Three distributions are shown for each component: all 60 cycles, the 10 fastest-burning cycles, and the 10 slowest-burning cycles.

In the text

Figure 15 Distributions of nine components of the velocity gradient tensor in the volume centered at the spark plug at SOIgn (342 °aTDCe) for two different LES models. The distributions are sampled at all volume nodes in the region for the 60 simulated cycles (total sample size of 370k). Three distributions are shown for each component: all 60 cycles, the 10 fastest-burning cycles, and the 10 slowest-burning cycles. This is the same as Figure 14, but plotted using a logarithmic scale.

In the text