A Diagnostic Approach to Assess the Effect of Temperature Stratification on the Combustion Modes of Gasoline Surrogates

ABSTRACT

Thermal stratification may switch the combustion mode from deflagration to spontaneous (auto-ignition) in spark-ignition engines leading to knock. Despite the available numerical and experimental works, an analytical and systematic method on the role of local thermal stratification on the combustion mode is still missing. Particularly, the effects of heat diffusion before ignition in negative temperature coefficient (NTC) chemistry of gasoline surrogates are typically ignored. In this study, an a-priori diagnostics tool is provided to separate the deflagration and auto-ignition combustion modes by considering the diffusion effects, based on two parameters: stratification wavelength ($\lambda$) and amplitude ($\delta$). The diagnostics tool is an extension of Zeldovich’s theory to transient problems by solving the diffusion equation and considering the flame and ignition timescales. It is found that 1) The theory is valid against one-dimensional numerical simulations under different average temperatures and pressures. 2) NTC chemistry promotes spontaneous ignition at both low and high pressures. 3) In the presence of NTC chemistry, the transition region between the two combustion modes is broadened. A third blended mode is observed (spontaneous ignition assisted flame) with front speeds $\approx 6S_l$. 4) Finally, estimation of the knock propensity from the generated maps is related to the surrogates’ octane sensitivities.

KEYWORDS:

• Combustion modes
• temperature stratification
• gasoline surrogates
• spark ignition engines
• knock

Introduction

The pathway toward using environmentally friendly fuels and exploiting higher thermal efficiencies in internal combustion engines has led to developments in fuel composition and combustion strategies in both Spark Ignition (SI) and Compression Ignition (CI) engines. Concerning the SI (gasoline fueled) engines, the major parameters affecting the engine efficiency are engine knock and emissions (Sarathy, Farooq, Kalghatgi 2018; Wang, Liu, Reitz 2017). While deflagration is the preferred mode of combustion to avoid knock in SI engines, engine knock is generally induced by auto-ignition of the end-gas and the subsequent pressure wave generation (Wang, Liu, Reitz 2017). Therefore, it is of paramount importance to identify between the two modes of combustion in SI engines: deflagration and auto-ignition.

In general, combustion initiated from a hot spot can be categorized into five modes; thermal explosion, supersonic auto-ignitive deflagration, developing detonation, subsonic auto ignitive deflagration and laminar deflagration (Gu, Emerson, Bradley 2003). According to the literature, the knocking process in SI engines can be attributed to the occurrence of the subsonic auto-ignitive deflagration mode (Bates et al. 2016; Wei et al. 2021). The aforementioned combustion mode is associated with slight or moderate knock in SI engines (Bates et al. 2016). Moreover, the intermittent generation of pressure waves and the induced ignition characterized by deflagration to detonation transition can result in severe knock and, in extreme cases, super knock (Wang et al. 2015; Wei et al. 2021; Zhou et al. 2020).

Combustion modes in SI engines may be affected by the fuel composition and thermal stratification. For instance, while bio-derived gasoline (bio-gasoline) is reported to reduce CO${ }_2$ emissions, it can enhance knock resistance of the fuel mixture (Cheng et al. 2020). On the other hand, thermal stratification effects on the combustion mode have been previously studied by different researchers as reviewed in the remainder of this section. In this work, a systematic approach is proposed to identify between two modes of combustion, i.e. deflagration and auto-ignition, with various thermal stratification levels using a simple analytical method, originally developed by Karimkashi et al. (2020). To validate the analytical model, hundreds of one-dimensional (1D) numerical simulations under different stratification levels and for two different gasoline surrogates are conducted.

The simplest and most commonly used gasoline surrogates are Primary Reference Fuels (PRF), which are binary mixtures of $n$-heptane and iso-octane. These mixtures are used as a reference for measuring the fuel anti-knock quality expressed as Research Octane Number (RON) and Motor Octane Number (MON) test conditions in the Cooperative Fuel Research (CFR) engine (ASTM 2013, 2014). With relevance to this work, a binary mixture of iso-octane/n-heptane (87%/13% volume, PRF 87) and a ternary mixture of n-heptane/iso-octane/ethanol (18% /62% /20% volume, PRF-E) are considered. Effects of chemical composition of the gasoline surrogates on their combustion characteristics are numerically addressed in the literature by correlating their ignition delay time (IDT) to the knock index. For instance, Naser et al. (2017) correlated the Octane Index (OI) and RON and MON numbers of various PRFs and Toluene Primary Reference Fuels (TPRFs) to their IDT at 21.3 bar and 770 and 850 K. Mehl et al. (2011b) correlated the anti-knock index of gasoline surrogates to their IDT at 25 bar and 835 K. Badra et al. (2015) expanded the correlation for TPRFs octane numbers using IDTs at 850 K and 50 bar. However, a key difference between the measured fuel knock quality and the calculated IDTs is that zero-dimensional (0D) simulations assume a premixed, fully homogeneous gas phase reactions.

Combustion mode classification of a reactivity stratified mixture is shown to be a beneficial approach for studying the transient flame dynamics under engine relevant conditions. For instance, the influence of thermal stratification on the auto-ignition of hydrogen/air mixtures under Homogeneous Charge Compression Ignition (HCCI) conditions is addressed using Direct Numerical Simulation (DNS) (Chen et al. 2006; Hawkes et al. 2006; Sankaran et al. 2005). Sankaran et al. (2005) derived a predictive criterion for deflagration and spontaneous ignition modes based on the ratio of the reaction front propagation speeds. The criterion has been found to be consistent with the results from DNS in Refs. (Chen et al. 2006; Sankaran et al. 2005). Chen et al. (2006) developed diagnostic measures for analyzing the mode of combustion. They have found that the prevalence of deflagrative versus spontaneous ignition modes depends on the local thermal stratification. Hawkes et al. (2006) further examined the effect of initial amplitude and wavelength of temperature fluctuations and turbulence intensity on ignition front propagation. Their DNS results demonstrated that higher amplitudes of temperature fluctuations can increase the prevalence of deflagration.

The effect of thermal stratification, turbulence, and spark ignition on premixed iso-octane/air and $n$-heptane/air mixtures under HCCI and Spark-Assisted Compression Ignition (SACI) conditions was studied by Yoo et al. (2013, 2011). DNS along with reduced mechanisms were utilized on a two-dimensional (2D) computational domain while the average temperature was kept higher than the Negative Temperature Coefficient (NTC) range. It was demonstrated that high levels of temperature stratification induce a deflagration mode in reaction fronts and thus result in smooth mean heat release rate. Dai et al. (2015) performed a one-dimensional (1D) numerical study on reaction front propagation modes of $n$-heptane/air mixture at initial temperatures within and below the NTC range inside a hotspot. Their results indicate that cool spots may also generate knock in engines when the engine is operated close to the fuel NTC range. Recently, Nogawa and Terashima (2020) studied the effect of global temperature stratification on the ignition modes in a constant volume reactor. In their study, $n$-heptane and $n$-butane/air mixture were chosen as NTC and non-NTC fuels. They observed that an increase in initial pressure and occurrence of NTC can both increase the possibility of heavy knocking.

Moreover, 1D study of Gupta, Im, and Valorani (2013) on $n$-heptane/air mixture showed that under relatively large temperature stratification levels and at non-NTC conditions, there could be a transition from deflagration to spontaneous ignition. Further, the mixture at the NTC range obtained two-stage ignition while the first ignition was mostly at deflagration regime and the secondary ignition was almost always at the spontaneous ignition mode.

With relevance to the fuel mixtures investigated in the present study, Luong et al. (2013) studied the effects of thermal stratification on the ignition of different PRF/air mixtures under HCCI conditions. The 2D DNS results revealed that at high levels of thermal stratification, the effect of PRF composition vanishes by observing similar deflagration reaction front propagation between the different mixtures. Recently, ignition characteristics of a mixture of PRF/air/exhaust gas recirculation (EGR) under reactivity controlled compression ignition (RCCI) and direct dual-fuel stratification (DDFS) was studied by Luong et al. (2017). The 2D DNS results showed that DDFS combustion results in longer combustion duration and lower peak of Heat Release Rate (HRR) compared to the RCCI combustion. More recently, Karimkashi et al. (2020) presented a theoretical relationship for predicting spontaneous and deflagrative reaction front propagation for locally stratified dual-fuel mixtures under RCCI engine-relevant conditions. For their test case of methanol/n-dodecane blends, the theoretical analysis was able to predict the phase border between the combustion modes observed in 1D chemical kinetics simulations. Furthermore, they developed the theory for turbulent cases and examined the effect of convection on the combustion modes using fully-resolved 2D simulations.

According to the above-mentioned literature review, experimental and numerical studies on effects of thermal stratification on the combustion mode for fuels and surrogates utilized in SI engines are quite extensive. However, an analytical model to explore thermal stratification effects on the combustion mode in a systematic way is still lacking. On the other hand, numerical approaches for evaluating the combustion characteristics of gasoline surrogates (e.g. knock tendency) mainly concern 0D kinetics modeling. In particular, the effect of diffusion in homogenizing the temperature prior to ignition is not addressed in the literature due to the choice of length scales such that the thermal diffusion is not active. Therefore, a systematic and comparative study on the effect of thermal stratification levels on the combustion mode should be provided for practical gasoline surrogates in SI engines with effects of diffusion included.

Therefore, in order to address the research gap outlined above and in continuation of the authors’ previous study (Karimkashi et al. 2020), here, we study a simplified 1D model problem with relevance to the unburnt fuel-oxidizer mixture. Next, we motivate the studied simulation setup. In a spark-ignition (SI) engine a premixed fuel-oxidizer mixture is ignited by a spark plug leading ideally to a flame propagation through a turbulent flow field. Such a setup is illustrated in Figure 1a) showing a schematic of an engine cylinder. One of the illustrated temperature cut-planes is taken from the datasets of our previous study by Masouleh et al. (2018), where we carried out large-eddy simulations (LES) of flame propagation in a SI model engine using the G-equation combustion model assuming a fully premixed fuel-oxidizer mixture. However, as noted, the temperature stratification can be strong near the cylinder walls. In such circumstances, laminar flame speed and ignition front speed can be strongly dependent on the local conditions, thereby affecting the local combustion progress. In the present study, we investigate different ignition scenarios of the end-gas using the model problem of sinusoidally varying temperature fluctuation as depicted schematically in Figure 1b). The model could effectively correspond to combustion process in a locally stratified region as depicted in Figure 1a).

Figure 1. Schematic presentation of a) SI engine cylinder and temperature distribution in the unburnt gas on a cut-plane inside the cylinder based on previous large-eddy simulations by the authors (Masouleh et al. 2018), b) model problem in the present study.

We note that the model problem in this study provides details on local combustion phenomena at relatively small length-scales using detailed chemistry and at high grid resolutions. This information can be utilized alongside CFD simulations where using detailed chemistry and high grid resolutions are inhibitive. The beta-curve can be also used as an a-priori tool to find the interesting conditions for CFD simulations.

This study aims to:

1. ) Extend the theoretical analysis presented for mixture stratification under compression ignition engines in Ref. (Karimkashi et al. 2020) to temperature stratification conditions in SI engines,

2. ) Analyze the effect of temperature stratification on the combustion mode of two different gasoline surrogates (PRF and PRF-ethanol) under normal and boosted engine relevant conditions, using a recently developed computational framework,

3. ) Explore the ignition front characteristics under different thermal stratification levels for the considered gasoline surrogates,

4. ) Provide a regime diagram in order to assess the combustion characteristics (knock tendency) of the gasoline surrogates, and

5. ) Formulate a diagnostic approach enabling the generation of predictive maps for various fuel surrogates under thermal stratification.

Methodology

Building upon the previous work by Karimkashi et al. (2020) on the combustion mode analysis of locally stratified mixtures, a similar strategy is taken to characterize the modes of combustion for different levels of initial temperature stratification. The present approach includes a theoretical analysis followed by hundreds of 1D numerical simulations for two different gasoline surrogates under SI engine-relevant conditions. The unsteady 1D simulations evaluate the validity of the theoretical analysis while providing details on the transient combustion characteristics. The details about these steps are discussed in the following.

Theoretical analysis

The theoretical approach developed by Karimkashi et al. (2020) is an a-priori tool to evaluate the prevelance of auto-ignitive versus deflagrative front propagation in locally stratified mixtures relevant to compression ignition engines. This theory, called $\beta$-curve analysis, is slightly modified herein to be utilized for the combustion mode analysis in homogeneous mixtures with different thermal stratification levels under SI engine-relevant conditions. In the following, the $\beta$-curve theory is briefly introduced and it is accordingly modified for thermal stratification.

Sankaran et al. (2005) defined a parameter, $\beta$, as an indicator for dominance of the spontaneous ignition versus deflagration regimes in a mixture with temperature stratification,

(1) $\beta =C_{\beta }\frac{S_l}{S_{ign}},$(1)

where $C_{\beta }$ is a coefficient of the order of unity, $S_l$ is the laminar flame propagation speed, and $S_{ign}$ is the front propagation speed after the initial ignition in the mixture. Therefore, $\beta >1$ represents the deflagration regime whereas $\beta <1$ implies the spontaneous ignition mode, and $\beta$=1 is the borderline between the two combustion modes. According to (Zeldovich 1980), $S_{ign}$ is equal to the inverse of IDT gradient for the premixed reactants. Assuming a homogeneous mixture, IDT is a function of temperature fluctuations only and $S_{ign}$ can be written as

(2) $S_{ign}=\frac{1}{|\mathrm{\nabla }\tau |}={\left(\frac{\mathrm{\partial }\tau }{\mathrm{\partial }T}\frac{\mathrm{\partial }T}{\mathrm{\partial }x}\right)}^{-1},$(2)

wherein $\tau$ represents IDT. Here, $\frac{\mathrm{\partial }\tau }{\mathrm{\partial }T}$, represents the ignition delay sensitivity to the temperature gradient. This term can also be interpreted as the slope of IDT when plotted against temperature ($\alpha$). The value of $\alpha$ changes depending on the fuel composition, average temperature, and pressure of the mixture. In Eq.(2(2) $S_{ign}=\frac{1}{|\mathrm{\nabla }\tau |}={\left(\frac{\mathrm{\partial }\tau }{\mathrm{\partial }T}\frac{\mathrm{\partial }T}{\mathrm{\partial }x}\right)}^{-1},$(2) ), $\frac{\mathrm{\partial }T}{\mathrm{\partial }x}$ is the temperature gradient with respect to the domain length, which reduces during the induction period (before auto-ignition) due to thermal diffusion.

Similar to the mixture stratification definition in Ref. (Karimkashi et al. 2020), the initial thermal stratification can be characterized using the spatial lengthscale ($\lambda$) and amplitude ($\delta$) in a presumed sinusoidal initial temperature stratification around an average temperature of $T_{ave}$, given by $T(x,t=0)=\delta sin(kx)$ in a 1D domain with periodic boundary conditions, i.e. $T(x=0,t)=T(x=\lambda ,t)$. Here, k =$\frac{2\pi }{\lambda }$ is the wavenumber. Karimkashi et al. (2020) provided an analytical solution for the 1D diffusion equation in ($\delta$, $\lambda$) space in order to estimate the diffusion effect with time evolution. Here, a similar approach is taken and the 1D heat diffusion equation provided below is solved by assuming a constant heat diffusivity and unity Prandtl number.

(3) $\frac{\mathrm{\partial }T}{\mathrm{\partial }t}=\nu \frac{{\mathrm{\partial }}^{2}T}{\mathrm{\partial }{x}^2}.$(3)

According to the provided details in Ref. (Karimkashi et al. 2020), Eq.3(3) $\frac{\mathrm{\partial }T}{\mathrm{\partial }t}=\nu \frac{{\mathrm{\partial }}^{2}T}{\mathrm{\partial }{x}^2}.$(3) is solved and the temperature gradient at the instant of ignition is written as

(4) ${\left.\frac{\mathrm{\partial }T}{\mathrm{\partial }x}\right|}_{t=\tau }=\delta kcos(kx)exp(-k^2\nu \tau ),$(4)

from which the maximum value of the temperature gradient in the domain can be written as

(5) ${\left({\left.\frac{\mathrm{\partial }T}{\mathrm{\partial }x}\right|}_{t=\tau }\right)}_{max}=\delta kexp(-k^2\nu \tau ).$(5)

Substituting Eq. 5(5) ${\left({\left.\frac{\mathrm{\partial }T}{\mathrm{\partial }x}\right|}_{t=\tau }\right)}_{max}=\delta kexp(-k^2\nu \tau ).$(5) into Eq. 2(2) $S_{ign}=\frac{1}{|\mathrm{\nabla }\tau |}={\left(\frac{\mathrm{\partial }\tau }{\mathrm{\partial }T}\frac{\mathrm{\partial }T}{\mathrm{\partial }x}\right)}^{-1},$(2) gives

(6) $S_{ign}{|}_{t=\tau }={\left(\alpha \delta kexp(-k^2\nu \tau )\right)}^{-1}.$(6)

Further, substituting Eq. 6(6) $S_{ign}{|}_{t=\tau }={\left(\alpha \delta kexp(-k^2\nu \tau )\right)}^{-1}.$(6) into Eq. 1(1) $\beta =C_{\beta }\frac{S_l}{S_{ign}},$(1) gives

(7) $\beta =C_{\beta }{S}_l\alpha \delta kexp(-k^2\nu \tau ).$(7)

Finally, considering $k=2\pi /\lambda$, the relation between $\lambda$ and $\delta$ can now be written as:

(8) $\delta =\frac{\beta \lambda exp(4{\pi }^2\nu \tau {\lambda }^2)}{C_{\beta }{S}_l\alpha 2\pi }.$(8)

Equation 8(8) $\delta =\frac{\beta \lambda exp(4{\pi }^2\nu \tau {\lambda }^2)}{C_{\beta }{S}_l\alpha 2\pi }.$(8) is the modified $\beta$-curve equation in Ref. (Karimkashi et al. 2020) for thermal stratification. Inserting $\beta =1$ into Eq.8(8) $\delta =\frac{\beta \lambda exp(4{\pi }^2\nu \tau {\lambda }^2)}{C_{\beta }{S}_l\alpha 2\pi }.$(8) provides the borderline between the deflagrative and spontaneous ignition modes in the ($\delta$,$\lambda$) space. The studied length scales are normalized by the diffusion length scale defined as ${\lambda }_{diff}=2\pi \sqrt{\nu \tau }$. A schematic plot of the $\beta$-curve is shown in Figure 2, wherein $\beta =1$ theoretically separates the two combustion modes from one another. Other $\beta$ values demonstrate the relative variation of the ignition front speed compared to the laminar flame speed, where $\beta >$ 1 ($<$ 1) represents dominance of deflagration (spontaneous) mode.

Figure 2. Schematic presentation of combustion regime diagram, the solid line ($\beta$=1) separates the spontaneous ignition and the deflagration modes. Different theoretical $\beta$ values (C${ }_{\beta }$=1) are also plotted (dashed lines). In practice, as noted also in the present work, the transition region around the $\beta$-curve may pose features from both combustion modes.

In the notation above, $\tau ,S_l$, and $\alpha$ are functions of $T_{ave}$, pressure, and the fuel composition. In the present study, values of ignition delay time ($\tau$) and laminar flame speed ($S_l$) are calculated using the 0D homogeneous reactor model and 1D premixed flame model of Cantera (Goodwin et al. 2018), respectively. The generated IDT vs temperature curves are used to calculate $\alpha$ at $T_{ave}$ as an estimation for local ignition delay sensitivity to temperature variations (Karimkashi et al. 2020). The transport properties library of Cantera has also been used for computation of the mixture viscosity.

Transient 1D simulations

The numerical framework utilized here is similar to the authors’ previous studies (Kahila et al. 2019; Karimkashi et al. 2020). An open source computational fluid dynamics (CFD) package OpenFOAM (Weller et al. 1998) is used to solve the 1D auto-ignition/flame propagation problem. The compressible Navier Stokes equations are coupled with species mass fraction transport equations along with the PIMPLE algorithm utilized for the pressure-velocity coupling. The chemical source terms are calculated using direct integration of finite rate chemistry. The customized chemistry solver, pyJac (Niemeyer, Curtis, Sung 2017), is used to evaluate the reaction rate coefficients. Additionally, pyJac offers an analytical Jacobian matrix formulation for solving the underlying chemical system of ODEs. For further details on the implementation and improvements over the standard OpenFOAM solvers, the reader is referred to (Kahila et al. 2019). It is noteworthy that the currently used customized solver has been shown to be able to capture the compression effects after the first ignition (Kahila et al. 2018).

A 1D domain with periodic boundary conditions is considered resembling the constant volume combustion in engines. The size of computational domain in all simulations is equal to $\lambda$. Initialization of the temperature fluctuation is set within a range of specific amplitudes (0 $<$ $\delta$ $<$ 200) and length scales (0 $<$ $\lambda$ $<$ 2${\lambda }_{diff}$) together with a homogeneous fuel mass fractions distribution. As shown by Masouleh et al. (2018), in the highly simplified model engine, the local temperature stratification can be as high as $\delta \approx 50-200$ K close to the walls while becoming more uniform in the center of the cylinder.

A computational grid with uniform grid spacing of 0.05 $\times$ flame thickness (${\delta }_f$) is utilized to resolve the flame scales. The flame thickness is calculated using the data from the laminar flame simulations, ${\delta }_f=(T_b-T_u)/{\left(\frac{\delta T}{\delta x}\right)}_{max}$ with $T_b$ and $T_u$ being the product and reactant temperatures, respectively. In 1D simulations, IDT has been considered as the time when the maximum temperature exceeds the threshold value of 1500 K ($T_{max}$ $>$ 1500 K). In the transient OpenFOAM simulations, the time step after the initial ignition is dynamically adjusted so that CFL ${ }_{max}$ =0.1 constraint.

A schematic picture of the anticipated combustion phases is shown in Figure 3. According to the sinusoidal temperature distribution, the mixture is classified as two distinct zones with higher and lower reactivity (step I). The key process occurring before initial ignition (step II) is heat diffusion. Diffusion reduces the temperature heterogeneity from the initial time until the ignition time (the dashed lines in Figure 3, Step I). Step II is attributed to the instance of time when the zone with higher reactivity starts to ignite. It should be noted that during step II, the reactivity of the lower reactive part of the mixture increases compared to its initial state due to the low temperature chemistry (LTC) and heating compression of the initial ignition kernel. Two distinct scenarios are expected to happen after the initial ignition. In the first scenario, step III-a, auto-ignition happens in the lower reactive part of the mixture independently and slightly after the initial ignition. This will form a secondary ignition front and it consumes rest of the unburnt reactants, resembling a volumetric reaction process. The main reason for the occurrence of step III-a is that the lower and higher reactive sides have similar reactivity stratification at the instant of primary ignition. The reason for this event can be either a very high thermal diffusion rate prior to the first ignition or initially low levels of reactivity stratification. The initial reactivity stratification level is related to the initial temperature profile. On the other hand, occurrence of NTC also modifies the reactivity distribution. The NTC behavior of the fuel can result in equivalent IDTs for different temperature levels. Therefore, the fuel type (NTC or no NTC), average temperature, and pressure may affect the reactivity stratification levels. The second scenario (step III-b) is that the initial ignition kernel forms a flame front witch propagates through the lower reactive mixture side. Considering the periodic boundary conditions in this study, it appears that the two fronts are approaching each other in the lower reactive side of the mixture (Figure 3, step III-b). The incidence of such process is due to the high reactivity stratification along the domain at the instant of primary ignition. Therefore, similar to the steady state flame propagation, the diffusion of heat and active radicals through the flame front induces reactions in the lower reactivity side of the mixture.

Figure 3. Schematic presentation of combustion process, step I thermal diffusion, step II primary ignition, step III-a in the spontaneous ignition mode an ignition front propagates toward the unburnt mixture, and step III-b in the deflagration mode a premixed flame propagates toward the unburnt mixture.

Fuel reactivity distribution vs temperature

In the present study, two different gasoline surrogates are examined for the combustion mode analysis. A binary mixture of iso-octane/$n$-heptane (87%/13% volume, PRF 87) and a ternary mixture of $n$-heptane/iso-octane/ethanol (18% /62% /20% volume, PRF-E) are considered. PRF mixture has similar MON value to gasoline, while PRF-E mixture composition mimics an oxygenated gasoline with MON = 89 and RON = 95.1. The choice of PRF-E mixture enables taking into account the octane sensitivity observed in the commercial gasoline and then also helps to observe the improvement of OI by using oxygenated additives. Addition of ethanol decreases the NTC behavior of the mixture, yet it increases the ignition propensity in intermediate temperature regions. On the contrary, the presence of ethanol increases the resistance to ignition in low temperatures where knock typically occurs. However, at high temperatures, the ignition promoting or inhibiting effect of ethanol vanishes, and biogasoline ignites similarly to gasoline (Frassoldati et al. 2010).

In order to choose an appropriate chemical mechanism, mechanism validation is shown in Appendix A. Additionally, comparisons between various mechanisms is also shown in Appendix A. High pressure 0D constant volume homogeneous reactor and 1D premixed flame models in Cantera are used for the validation. Consequently, a 115 species skeletal TPRF mechanism developed by Stagni et al. (2016) for PRF simulations and a skeletal 171 species bio-gasoline mechanism (Ranzi et al. 2014) for the PRF-E simulations are selected.

We note that recent engine development has led to higher pressures at fixed temperatures, causing engines to work under RON-like conditions or even beyond that. The reciprocal relation of the fuel sensitivity and knock tendency of the fuels at high loads motivated investigation of PRF and PRF-E combustion modes at two different pressures, i.e. 20 bar and 50 bar. Figures 4a and 4b depict IDT versus temperature distribution of the investigated fuels at 20 and 50 bars, respectively. While both of the IDT profiles exhibit NTC at 20 bar pressure, there is mild NTC behavior for PRF-E at 50 bar. At 20 bar, Figure 4a, two average temperatures are considered for further examinations. Point I at $T_{ave}$=1000 K is in the high temperature region and point II at $T_{ave}$=830 K, is located at the NTC region for both fuels. At 50 bar, three different $T_{ave}$ at 1000 K (III), 900 (IV) and 830 (V) are selected for high, intermediate, and low temperature zones, respectively. In the following 1D simulations, these points are the average temperatures ($T_{ave}$) in cases I–V around which different temperature amplitudes ($\delta$) with various lengthscales ($\lambda$) are considered. Detailed information on the reference values calculated at $T_{ave}$ for these cases at points I–V are reported in Table 1. It is worth mentioning that according to the high reactivity of the mixture at the examined pressures and temperatures, the presented reference values for the flame speed and flame thickness are estimated values using Cantera (Krisman, Hawkes, Chen 2018).

Figure 4. PRF and PRF-E ignition delay time versus temperature profiles at P = 20 and 50 bar, black circles are marking the locations of the studied cases.

Table 1. Conditions for the five cases: pressure, average temperature, ignition delay time, flame speed, product temperature, and flame width.

	Fuel	P [bar]	T${ }_{ave}$ [K]	IDT[ms]	S${ }_l$ [m/s]	T ${ }_b$ [K]	${\delta }_f$ [mm]	${\lambda }_{diff}$ [mm]
Case I	PRF	20	1000	2.9	2.22	2700	0.015	0.8
	PRF-E	20	1000	2.2	2.1	2700	0.016	0.7
Case II	PRF	20	830	6.7	1.31	2600	0.019	1
	PRF-E	20	830	11	1.24	2600	0.02	1.3
Case III	PRF	50	1000	1.2	1.73	2700	0.007	0.32
	PRF-E	50	1000	0.9	1.63	2600	0.008	0.28
Case IV	PRF	50	900	1.6	1.3	2700	0.008	0.35
	PRF-E	50	900	2.2	1.2	2700	0.009	0.4
Case V	PRF	50	830	1.5	1.08	2600	0.009	0.31
	PRF-E	50	830	2.3	0.95	2600	0.01	0.39

Results and discussions

Combustion mode analysis

In this section, combustion mode analysis for different case studies defined in Table 1 are provided using 1D numerical simulations as well as the $\beta$-curve theory. As it was explained in section 2.1, $\beta >1$ and $\beta <1$ designate spontaneous ignition and deflagration combustion modes respectively. Moreover, $\beta \approx$ 1 (Eq 1)(1) $\beta =C_{\beta }\frac{S_l}{S_{ign}},$(1) describes a reaction zone where spontaneous ignition and deflagration have equivalent role in front propagation and $C_{\beta }$ is a coefficient of the order of unity. Moreover, it has been shown that Eq 8(8) $\delta =\frac{\beta \lambda exp(4{\pi }^2\nu \tau {\lambda }^2)}{C_{\beta }{S}_l\alpha 2\pi }.$(8) draws a theoretical borderline between the ignition and deflagration zones in ($\lambda$, $\delta$) space.

The following figures (Figures 5–9) depict the output of 1D simulations as well as the $\beta$-curves calculated using the average temperature and pressure with different values of $0.4\le C_{\beta }\le 3$. Each contour contains about $\approx$ 100 simulation points, in which $\lambda$ and $\delta$ are varied to achieve different initial sinusoidal temperature distributions. A temperature threshold (1500 K) is considered in order to distinguish the start of ignition for each case. After the initial ignition, a characteristic flame formation time is calculated for each case ($d{t}_f=\frac{{\delta }_f}{S_l}$) and it is added to the respective ignition time. Eventually, the contours have been created for the maximum temperature difference in the domain ($\mathrm{\Delta }{T}_{max}$) achieved for each case after the mentioned time period, i.e. at IDT + $d{t}_f$. The color scale in the contours is based on how close $\mathrm{\Delta }{T}_{max}$ is to a cutoff value of temperature difference, $\mathrm{\Delta }{T}_{th}=1200$ K, considered in all cases. The choice of $\mathrm{\Delta }{T}_{th}$ value, is based on the difference between the average temperature of the unburnt mixture (830 K $<$ T${ }_u<$ 1000 K) and burnt products (T${ }_b\approx$ 2600 K) of a steady state premixed flame, according to Table 1. The local temperature gradient is used to distinguish between a flame and an ignition front (Karimkashi et al. 2020; Sankaran et al. 2005). The regions where $\mathrm{\Delta }{T}_{max}\ll \mathrm{\Delta }{T}_{th}$ are considered to be governed by spontaneous ignition while those with $\mathrm{\Delta }{T}_{max}\ge \mathrm{\Delta }{T}_{th}$ are interpreted to be characterized by the deflagration mode. In the following, each figure corresponding to the regime diagrams of each case study for the two considered gasoline surrogates is analyzed, comparatively.

Figure 5. Combustion mode distribution of case I ($P$ = 20 bar and $T_{ave}$ = 1000 K), black and white regions represent spontaneous ignition and deflagration modes of a) PRF and b) PRF-E mixtures respectively. Colored lines display the $\beta$-curves calculated using different $C_{\beta }$ values, Eq.8(8) $\delta =\frac{\beta \lambda exp(4{\pi }^2\nu \tau {\lambda }^2)}{C_{\beta }{S}_l\alpha 2\pi }.$(8) . c) IDT distribution at the initial time versus T of PRF and PRF-E mixtures.

Figure 6. Combustion mode distribution of case II ($P$ = 20 bar and $T_{ave}$ = 830 K), black and white regions represent spontaneous ignition and deflagration modes of a) PRF and b) PRF-E mixtures respectively. Colored lines display the $\beta$-curves calculated using different $C_{\beta }$ values, Eq.8(8) $\delta =\frac{\beta \lambda exp(4{\pi }^2\nu \tau {\lambda }^2)}{C_{\beta }{S}_l\alpha 2\pi }.$(8) . c) IDT distribution at the initial time versus T of PRF and PRF-E mixtures.

Figure 7. Combustion mode distribution of case III ($P$ = 50 bar and $T_{ave}$ = 1000 K), black and white regions represent spontaneous ignition and deflagration modes of a) PRF and b) PRF-E mixtures respectively. Colored lines display the $\beta$-curves calculated using different $C_{\beta }$ values, Eq.8(8) $\delta =\frac{\beta \lambda exp(4{\pi }^2\nu \tau {\lambda }^2)}{C_{\beta }{S}_l\alpha 2\pi }.$(8) . c) IDT distribution at the initial time versus T of PRF and PRF-E mixtures.

Figure 8. Combustion mode distribution of case IV ($P$ = 50 bar and $T_{ave}$ = 900 K), black and white regions represent spontaneous ignition and deflagration modes of a) PRF and b) PRF-E mixtures respectively. Colored lines display the $\beta$-curves calculated using different $C_{\beta }$ values, Eq.8(8) $\delta =\frac{\beta \lambda exp(4{\pi }^2\nu \tau {\lambda }^2)}{C_{\beta }{S}_l\alpha 2\pi }.$(8) . c) IDT distribution at the initial time versus T of PRF and PRF-E mixtures.

Figure 9. Combustion mode distribution of case V ($P$ = 50 bar and $T_{ave}$ = 830 K), black and white regions represent spontaneous ignition and deflagration modes of a) PRF and b) PRF-E mixtures respectively. Colored lines display the $\beta$-curves calculated using different $C_{\beta }$ values, Eq.8(8) $\delta =\frac{\beta \lambda exp(4{\pi }^2\nu \tau {\lambda }^2)}{C_{\beta }{S}_l\alpha 2\pi }.$(8) . c) IDT distribution at the initial time versus T of PRF and PRF-E mixtures.

Regime diagram: case I

Figure 5 presents 1D simulation results for $T_{ave}$ = 1000 K and P = 20 bar. It is observed from the combustion mode diagrams, Figures 5a and 5b, that fluctuations with higher $\lambda$ and $\delta$ values promote deflagration mode compared to spontaneous ignition. The overall expansion of the deflagrative mode under higher levels of temperature stratification is also observed in previous studies using iso-octane/air and PRF/air mixtures (Luong et al. 2013; Yoo et al. 2013).

A particular characteristic phenomenon in the proximity of diffusion length scales is that increasing $\lambda$ results in less homogeneity due to the weakened diffusion process prior to the initial ignition and accordingly, higher reactivity stratification levels at the onset of ignition. However, above the considered scales, an increase in $\lambda$ reduces the temperature gradient and hence results in a more homogeneous mixture and stronger ignition.

It is worth noting that in this study we always keep the maximum wave lengths close to the value of the diffusion length scale ($\frac{\lambda }{{\lambda }_{diff}}\le 2.0$) to prevent the formation of intense pressure waves, and thereby initiation of the detonation phenomenon, at the instance of ignition in the domain.

Figures 5a and 5b depict that combustion mode distribution does not change along the horizontal axis when the amplitude increases to values higher than 100 K. This is due to the large difference between the IDT values of the minimum and the maximum temperature points (e.g. $\delta$ = 100, 730 K $<$ T $<$ 930 K), taking into account the monotonic shape of the IDT vs temperature profiles. The reactivity distribution of PRF and PRF-E versus $\delta$ at $T_{ave}$ = 1000 K is displayed in Figure 5c.

When comparing Figures 5a and 5b, both fuels are seen to behave similarly in terms of combustion mode distribution. This behavior may also be deduced from Figure 5c. The present result is inline with kinetic study of Frassoldati et al. (2010), in which they depicted the similar IDTs of gasoline and bio-gasoline surrogates under high temperature conditions.

Regime diagram: case II

Figure 6 presents 1D simulation results for $T_{ave}$ = 830 K and P = 20 bar conditions. It is observed that the PRF mixture is more prone to spontaneous ignition at similar stratification levels, which is directly related to the reactivity stratification differences between the two mixtures. According to Figure 6c, although the IDT profiles of both fuels show a considerable NTC behavior, the PRF-E mixture shows lower degrees of NTC compared to the PRF mixture. Therefore, reactivity stratification of the PRF mixture is considerably smaller than that of the PRF-E mixture at the considered temperature and pressure range. Accordingly, the extent of spontaneous ignition is lower for the PRF-E compared to the PRF.

We note that according to Figure 4a, when approaching lower temperatures, at T${ }_{ave}$ = 830 K and P = 20 bar, the PRF mixture tends to be more reactive than the PRF-E mixture. However, conducted research shows that correlating IDT values and knock indexes, e.g. RON, for gasoline and bio-gasoline surrogates can result into significant errors (Naser et al. 2017; Singh et al. 2019). Although longer IDTs of the PRF-E at low temperatures can decrease the ignition propensity of the mixture (Cheng et al. 2020), it allows diffusion to neutralize the initial temperature stratification (larger ${\lambda }_{diff}$ values) and therefore promote spontaneous ignition. By analyzing the events after initial ignition and using the dimensionless length-scales ($\frac{\lambda }{{\lambda }_{diff}}$), the effect of the ignition delay time (IDT) offset on the results has been normalized for all cases. Therefore, as it is apparent from Figures 6a and Figure 6b, at high values of $\delta$, both fuels have the same combustion mode distribution.

A comparison between Figures 6a and Figure 5a, shows that spontaneous ignition region is advanced toward the higher $\delta$ and $\lambda$, implying a higher spontaneous ignition propensity for the PRF. At $T_{ave}$ = 830 K, the NTC behavior reduces the variation of reactivity versus amplitudes and results in higher chance of spontaneous ignition, within the range shown in Figure 6c. Considering the differences between Figures 6b and Figure 5b, a similar conclusion can be made for the PRF-E fuel. The reduced reactivity variation results in higher significance of spontaneous ignition mode for the PRF-E fuel.

In addition, it is observed that the theoretical $\beta$-curve region between the spontaneous ignition and the deflagration modes is wider compared to the case I, especially for the PRF-E mixture. This may be connected to a situation where both combustion modes are significant (Bates et al. 2016). Here, when both modes are important, it is labeled as spontaneous ignition assisted deflagration. The topic is analyzed more in sections 3.1.6 and 3.2.

Regime diagram: case III

Figure 7 presents 1D simulation results for $T_{ave}$ = 1000 K and P = 50 bar. Although the combustion mode profiles look similar to Case I, according to Table 1, ${\lambda }_{diff}$ obtains lower values for both fuels. In fact, considering lower stratification levels (upper branch) the deflagration zone is now limited to smaller length scales compared to Case I. With the pressure increase to 50 bar and at high average temperature (T${ }_{ave}$ = 1000 K), Figure 7c shows that the PRF IDT profile exhibits a mild NTC at 900 K $<$ T $<$ 950 K. The displacement of the NTC region toward higher temperatures with increasing pressure is related to the competition between $H_2{O}_2$ decomposition and low-temperature chain routes (Lee et al. 2017). Compared to the case I (Figure 5a), the spontaneous ignition zone shifts to higher $\lambda$ and $\delta$ values due to the appearance of the NTC in the PRF profile. The PRF and PRF-E combustion mode distributions, Figures 7a and 7b, can be compared with a similar analogy to case I. Increasing $\lambda$ and $\delta$ results in higher reactivity stratification levels at the onset of ignition which leads to enhanced deflagration mode for both mixtures. However, at comparatively lower $\lambda$ and $\delta$ values, the spontaneous ignition mode is more dominant for the PRF due to the presence of NTC chemistry.

Moreover, it is observed from Figures 7a and 7b that the combustion mode distribution does not change by increasing $\delta$ to higher values than 150 K. This increase in the maximum effective amplitude compared to case I is a result of the lower gradient of reactivity versus temperature, Figure 7c, under elevated pressures.

Regime diagram: case IV

Figure 8 presents 1D simulation results for $T_{ave}$ = 900 K and P = 50 bar. Figure 8c shows that by changing the average temperature to 900 K, the PRF mixture exhibits a considerable NTC behavior leading to lower reactivities at higher temperatures. However, the PRF-E profile shows a mild NTC behavior. As a result, the combustion mode distribution of the PRF-E mixture remains fairly similar to the case III, while the spontaneous ignition zone for the PRF mixture proceeds toward the higher amplitudes. The observed combustion mode distribution for the two surrogates is in line with the fact that the PRF (zero sensitivity) is more prone to knock than a sensitive (real gasoline or ethanol blended PRF) fuel at higher pressures (Szybist and Splitter 2017).

Regime diagram: case V

Figure 9 presents 1D simulation results for $T_{ave}$ = 830 K and P = 50 bar. Shifting to lower average temperatures (T ${ }_{ave}$ = 830 K), the amount of NTC decreases for the PRF mixture, c.f. Figure 9c. Thereby, spontaneous ignition zone recedes to lower amplitudes compared to case IV, although it is still larger than in case III. However, as it is observed in Figure 9b, there is no noticeable change in the combustion mode distribution of the PRF-E.

The fact that combustion mode distribution of the PRF-E does not change significantly with changing the average temperature from 1000 K to 830 K under elevated pressures (cases III–V), illustrates the better performance of the ethanol blend in downsized-boosted engines with lower compression temperatures and higher pressures than the RON test. In addition, it is again observed that the theoretical $\beta$-curve region between the spontaneous ignition and the deflagration modes is wider compared to e.g. case I. This is connected to the previously mentioned situation where both combustion modes have significance (Bates et al. 2016), i.e. spontaneous ignition assisted deflagration.

Estimation of knock propensity

As mentioned earlier, knock tendency of the fuels is attributed to their ignition properties at the end gas. Although we do not simulate actual knock in the present study, the present numerical results still offer a possibility to estimate certain metrics to assess the mixture knock propensity. Here, we estimate the knock propensity of the fuels at each studied condition by calculating the relative area of the spontaneous ignition zone from 1D simulations. The area is then compared to the entire area of the $\lambda$-$\delta$ regime diagrams as percentages, c.f. Figures 5–9a and b. Table 2 summarizes the estimated knock propensity and NTC levels in cases I–V.

Table 2. NTC levels and estimated knock propensity of cases I–V based on the percentage of spontaneous ignition zone from 1D simulations. Additionally, the table gathers certain findings from cases I–V.

	Fuel	P[bar]	T${ }_{ave}$ [K]	NTC level	Knock Propensity [%]	Remarks
Case I	PRF	20	1000	No	23	narrow transition region
	PRF-E	20	1000	No	24
Case II	PRF	20	830	High	47.6	broad transition region
	PRF-E	20	830	High	33.3
Case III	PRF	50	1000	Mild	33	narrow transition region
	PRF-E	50	1000	No	23.7
Case IV	PRF	50	900	High	56.5	broad transition region
	PRF-E	50	900	Mild	30.4
Case V	PRF	50	830	Mild	49.8	broad transition region
	PRF-E	50	830	No	33.8

Similar knock tendency of the fuels under Case I conditions is also inline with their similar MON numbers, while the higher knock resistance of the PRF-E at high pressures and low temperatures (Cases IV and V) conforms with its higher RON number. Therefore, the enhancing effect of fuel sensitivity on knock resistance can be observed by comparing the estimated knock propensity of the PRF and PRF-E mixtures in Cases I and IV–V.

From Table 2 it is observable that the knocking tendencies are correlated with the NTC levels of the fuels. However, it is also worth noting that under certain levels of stratification (e.g. at high $\delta$ or $\lambda$ values) the knock propensities can be different from the respective estimation in Table 2. This quantification is not possible to be extracted from chemical kinetic analysis alone. Here the $\beta$-curve analysis enables taking the diffusion effect on temporal and spatial variation of temperature into account.

Discussion on the $\beta$-curve results

Figures 5–9a and b show the $\beta$-curves plotted using 0.4 $<C_{\beta }<$ 3.0. The idea of using different $C_{\beta }$ values is to identify the regions where ignition and deflagration modes have a similar significance. Therefore, the width of the confined zone between the upper and lower limits shows the rate of transition between the combustion modes.

For instance, the results for high temperature cases, I and III, show that the sensitivity of the $\beta$-curve location to $C_{\beta }$ is marginal. This implies that the transition from spontaneous ignition to deflagration combustion mode happens in a rather narrow region. This observation is in accordance with the similar IDT profiles and the clear trend between the temperature and reactivity of the mixture in cases I and III. However, $\beta$-curves for cases II, IV and V show high sensitivity to the values of $C_{\beta }$. This behavior implies that there is a broader range where both ignition and deflagration modes have a comparable significance (Bates et al. 2016). In this case the two modes could affect each other, i.e. increased propagation speed due to auto-ignition (spontaneous ignition assisted deflagration). Furthermore, when both modes are important, the rate of transition between the combustion modes is slower in low temperature conditions than in the high temperature cases. The current behavior is also in accordance with the moderate change observed in IDT vs temperature distribution in these cases, which is correlated to low and intermediate temperature heat release.

We note that the comparison between the theoretical curves and contours in Figures 5–9a and b show that the defined borderlines using the $\beta$-curve theory conforms with the results from 1D numerical simulations.

To reproduce a similar diagnosis to this study, the following steps should be taken:

(1) Finding 0D IDT and steady-state 1D premixed flame data (S${ }_l$, $\delta$, $\tau$, $\alpha$,$\nu$) at specific T${ }_{ave}$ values.

(2) Evaluating ${\lambda }_{diff}$ values and then drawing $\beta =1$ curve from Eq.(8)(8) $\delta =\frac{\beta \lambda exp(4{\pi }^2\nu \tau {\lambda }^2)}{C_{\beta }{S}_l\alpha 2\pi }.$(8) in the ($\lambda$, $\delta$) space.

(3) Performing a sensitivity test for different values of $C_{\beta }$ in the order of unity.

(4) Determining the surface area below $\beta =1$ curve to find the prevalence of ignition versus deflagration modes.

(5) As a complementary step, 1D simulations can be conducted, with the same approach discussed in the present paper, with the points chosen in proximity of the $\beta =1$ curve in ($\lambda$, $\delta$) space.

Note that if T ${ }_{ave}$ lies in the NTC zone the value of $\alpha$ becomes sensitive to the choice of T${ }_{ave}$ and chemical mechanism. Therefore, while conducting step 5 is strongly advised, achieved results can only show a qualitative comparison of the two combustion modes.

We note that the presented 1D diagnosis approach, provided with information on the average temperature, pressure, temperature stratification levels and fuel surrogates is able to predict the conditions that should be avoided to reduce the knock tendency. In practice the temperature stratification levels need to be estimated based on its correlation to turbulence levels and engine geometry (Kalghatgi et al. 2009). To estimate turbulence levels, it has been shown by Masouleh et al. (2019) that turbulence levels are correlated with engine RPM in practice. Furthermore, temperature stratification level varies inside the cylinder in real engines, therefore analyzing the sensitivity of knocking tendency to the temperature stratification is required which is also considered in the present diagnostic approach. Moreover, a possible extension of the present numerical framework to 2D, similar to the work by Karimkashi et al. (2020), may provide additional information on the effect of turbulent diffusion length scales in the current problem. Finally, a full 3D simulation is capable of capturing the phenomenon such as turbulence, cycle to cycle variation and geometrical effects. However, such simulations normally are limited by computational resources.

Analysis of the transient flame behavior

The intent of this section is to investigate temporal evolution of the ignition front in two cases from the previous section with different combustion mode distributions. In the previous section, we observed that there is a noticeable difference between the combustion mode distribution of the PRF-E and PRF mixtures at 900 K average temperature and 50 bar conditions. Figure 10 indicates the selected points for the transient study in this section. Points with analogous names in Figures 10a and 10b have identical normalized length scales and amplitudes. Point 1 is located in the deflagration zone for both mixtures with different distances from the spontaneous ignition zone. Point 2 is located at different combustion regimes, spontaneous mode for the PRF mixture and deflagration mode for the PRF-E mixture. Finally, point 3 is located in the deflagration zone for both mixtures at high values of $\delta$, similarly away from the spontaneous ignition zone. White lines in Figures 10a and 10b are associated with the theoretical $\beta$-curves for PRF and PRF-E mixtures with $C_{\beta }$ of 0.4 and 3, respectively. The $C_{\beta }$ s are chosen based on the conformity of the theoretical $\beta$-curves with 1D simulation results in Figure 8.

Figure 10. Locations of the selected points from (a) PRF and (b) PRF regime diagram for transient analysis at 900 K and 50 bar, white line represents theoretical $\beta$-curves calculated with $C_{\beta }$ values of 0.4 and 3, respectively.

We note that the estimated laminar flame speed for the current configuration ($T_{ave}$ = 900 K and $P$ = 50 bar) is 1.3 m/s, Table 1, for both the PRF and PRF-E mixtures. Therefore, the same front displacement speed is expected in the regions associated to the deflagration mode. Once the maximum temperature in the domain reaches 1500 K, it is considered as the initial ignition event with zero timing in the following discussion.

Figures 11–13 demonstrate transient profiles of the temperature and HRR for the selected points in Figure 10. In order to demonstrate the propagation speed, x-axis shows the real domain size and the profiles are plotted every 4 µs (except for point 2 for the PRF mixture). Thereby, using the domain length and the time period between the primary ignition (when $T_{max}>1500$ K) and the end of the process, the propagation speed can be estimated in each case.

Figure 11. Temporal evolution of temperature and HRR at point 1, $T_{ave}$ = 900 K and $P$ = 50 bar, demonstrating spontaneous ignition assisted deflagration for a) PRF and deflagration mode for b) PRF-E mixtures. Depicted profiles illustrate the propagation rates, initial and secondary ignitions.

Figure 12. Temporal evolution of temperature and HRR at point 2, $T_{ave}$ = 900 K and $P$ = 50 bar, demonstrating a) spontaneous ignition propagation for PRF mixture and b) spontaneous ignition assisted deflagration for the PRF-E mixture. Plotted profiles illustrate the occurrence of successive ignitions for the PRF mixture and the slower process for the PRF-E mixture. We note that propagation rates for both cases are higher than normal deflagration.

Figure 13. Temporal evolution of temperature and HRR at point 3, $T_{ave}$ = 900 K and $P$ = 50 bar, demonstrating deflagration mode for a) PRF b) PRF-E mixtures. Illustrated profiles show the propagation rates, initial ignition, and the consumption of reactants by approaching flames.

Transient flame analysis: point 1

As discussed earlier, the front propagation mechanism can either be spontaneous ignition or deflagration. In the former mechanism, homogeneously distributed reactions are driving a volumetric auto-ignition event, and the latter is a local and narrow reaction front. Therefore, one way to distinguish between these two modes is to study the heat release rate (HRR) evolution during the combustion process. In the spontaneous ignition mode, it is expected to observe almost evenly distributed HRR at different instances, whereas in the deflagrative mode, a wave like front at flame location is anticipated.

Figures 11a and 11b show transient temperature and HRR profiles at point 1 for the PRF and PRF-E, respectively. The HRR profiles show a flame-like structure with distinctive peaks indicating the front location. The initial time steps (t $<$ 8 $\mu s$) do not show a noticeable heat release in both of the PRF and PRF-E cases, which illustrates the initial ignition at early stages of flame formation. As the two fronts propagate toward each other, the remaining charge is heated by compression, resulting in the gradual HRR increase. In the final stages, the heated charge auto-ignites again (secondary ignition), resulting to decrease of the HRR and a sudden rise in the temperature profile.

Comparison between the front displacements during each time period indicates the difference in the front propagation mechanisms of the two fuels. From the shown profiles in Figures 11a and 11b, the front speed is estimated to be 9.67 m/s for the PRF mixture and 4.43 m/s for the PRF-E mixture. Based on these values, front propagation speeds are larger than the laminar flame speed in both mixtures. This higher propagation speed is expected due to the transient nature of the fronts compared to the steady state solution of the laminar flame. While for the PRF-E, front propagation speed is closer to the laminar flame speed, for the PRF mixture, this difference is larger. The high temperature of the unburned mixture ($T_{unburned}\ge 1000$ K) implies conditions without LTC chemistry, c.f. Figure 4, and the formation of spontaneous ignition assisted deflagration in the domain. The mentioned condition is defined as flame propagation enhanced by pre-ignition (HTC) of the reactants downstream of the flame (Krisman, Hawkes, Chen 2018). As a result, the front propagation can be classified as spontaneous ignition assisted deflagration and deflagration for the PRF and PRF-E mixtures, respectively. Considering the location of the point 1 in Figures 10a and 10b, the result also indicates that the closer the point is to the spontaneous ignition regime, the higher is the displacement speed.

Transient flame analysis: point 2

Figure 12 presents temporal distribution of temperature and HRR at point 2 for the PRF and PRF-E mixtures. Point 2 is located in the spontaneous ignition region for the PRF mixture and in the deflagration zone for the PRF-E fuel. Figure 12a shows that the temperature profile evolves in a rather homogeneous manner for the PRF mixture. As a result of the initial non-uniform temperature distribution, ignition starts at a higher temperature, forming a front at 4 $\mu s$. However, the unburnt part ignites instantaneously in the next time step. It is notable that the front formation is not necessarily an indicator of the deflagration mode. This is further discussed in the transport budget terms analysis section, where it is shown that the front observed here is an ignition front.

Transient profiles of the PRF-E mixture show similar trends with point 1 implying the occurrence of the deflagration front propagation. The estimated burning rate for the PRF and PRF-E are 120 m/s and 9.25 m/s, respectively. The former value is orders of magnitude higher than the laminar flame speed, which suggests the occurrence of spontaneous ignition while the latter value is in the same order of magnitude with laminar flame speed and implies an spontaneous ignition assisted deflagration mode.

Transient flame analysis: point 3

Figure 13 illustrates temporal evolution of the HRR and temperature at point 3. Transient profiles for both mixtures show a deflagration-like front propagation. Temperature profiles show that unlike point 1, the unburnt part of the mixture is consumed by deflagration mode, while it remained cold (T $<$ 1500 K) and hence secondary ignition (explained in Section 3.2.1) is avoided.

Moreover, estimation of the burning rate gives 3.94 m/s and 3.40 m/s front speeds for the PRF and PRF-E mixtures, respectively. The transient profiles and calculated reaction rates imply that the combustion modes are similar for both fuels at point 3 and it is classified as deflagration. The fact that the fuel reactivity effect on the combustion mode vanishes under high temperature stratification conditions has been also shown in a DNS study of lean PRF/air mixtures in Ref. (Luong et al. 2013).

The general trend observed in Figures 11–13 is inline with the results illustrated in Figure 10, confirming the validity of the assumptions made for the plotted regime diagrams. However, as noted before, front propagation is not solely representative of normal deflagration. This aspect has been investigated with transport budget analysis in the following section.

Transport budget terms analysis

Under the initial and boundary conditions considered in the present study, the convection term has a negligible contribution in the energy transport equation within the reaction zone of premixed flames (Gordon et al. 2007). Therefore, at normal deflagration mode, the reaction and diffusion terms are expected to be in balance, whereas for the spontaneous ignition mode reaction term dominates over the diffusion term (Karimkashi et al. 2020; Luong et al. 2017). Here, budget term analysis for the energy equation has been used to investigate combustion modes of the specified points in Figures 10a and 10b. Figure 14 compares the normalized diffusion and reaction terms of the points 1–3 for the PRF and PRF-E fuels at the same time instances inspected in the previous section. The diffusion and reaction terms in each figure are normalized with respect to the corresponding maximum value of the reaction term in each case.

Figure 14. Transport budget terms for points 1–3 of PRF and PRF-E mixtures, at $T_{ave}$ = 900 K and $P$ = 50 bar, representing the front propagation modes. Point 1: deflagration dominates (PRF-E) vs spontaneous ignition assisted deflagration (PRF), point 2: spontaneous ignition dominates (PRF) vs spontaneous ignition assisted deflagration (PRF-E), point 3: deflagration dominates (PRF and PRF-E).

A comparison between the budget terms at point 1 shows that for both fuels, diffusion and reaction terms have similar order of magnitude (reaction term has higher absolute values). Moreover, the coincidence of reaction and diffusion peaks suggests the existence of a thin reaction front which is a characteristic feature of the deflagration mode. Furthermore, the difference between the reaction and diffusion terms for the PRF mixture is higher than for the PRF-E mixture which explains the reason for the faster propagation of the front for the PRF mixture. This also implies a more significant deflagrative propagation for the PRF-E mixture. Moreover, this observation is in agreement with point 1 being farther away from the spontaneous ignition zone in the PRF-E mixture compared to the PRF in Figures 10a and 10b. Therefore, at point 1, the PRF mixture is classified as spontaneous ignition assisted deflagration, while the PRF-E is in the deflagration mode.

The budget terms at point 2 show two distinct combustion modes for the PRF and PRF-E mixtures. Using a similar analogy to the previous point, the combustion mode can be categorized as spontaneous ignition assisted deflagration for the PRF-E mixture. However, the PRF mixture profiles show a negligible contribution of the diffusion term while the reaction profiles are evenly distributed at most of the time instances. Furthermore, a closer look at the front shows that the diffusion term is two orders of magnitude smaller than the reaction term. This leads to the conclusion that the dominant combustion mode at point 2 for the PRF mixture is spontaneous ignition.

Comparison of the budget terms at point 3 shows that the contribution of the diffusion term increases in the energy transport with the flame evolution even though the reaction term has comparable values. Therefore, it can be concluded that the deflagration mode has a substantial role in the front propagation at point 3 for both mixtures.

It is noteworthy that, the budget term analysis at points 1 (for both mixtures) and 2 (only for the PRF-E mixture) illustrates the sequential occurrence of different combustion modes in the domain with time progress. Transition from the ignition dominated combustion mode (initial ignition) to the deflagration mode and again to the ignition mode (secondary ignition) can be deduced from the domination of the reaction term at initial and latest time steps. This trend is attributed to the initial ignition and secondary ignition due to the heated, compressed end-gas at the latest stages. The present investigation demonstrates that the occurrence of the secondary ignition depends on two factors. First, initial reactivity of the mixture, and second, the interaction between the approaching reaction fronts. The former factor can be controlled by changing the initial temperature stratification and fuel composition. However, the latter relies on the dynamics of flame propagation (Chen and Sohrab 1995).

Conclusions

The present work follows the numerical framework presented earlier by the authors on combustion mode analysis of a locally stratified dual-fuel mixture under CI engines relevant conditions. Here, the previous work is extended by studying the spontaneous ignition and deflagration front propagation in the end-gas using different characteristic length scales ($\lambda$) and amplitudes ($\delta$) of temperature stratification under SI engine relevant conditions. The work combined numerical and theoretical analyses. The theoretical results indicate a borderline between the spontaneous and deflagrative reaction front propagation regimes in $\lambda$-$\delta$ space by solving a one-dimensional, time dependent, diffusion-reaction problem with initially sinusoidal temperature stratification. Furthermore, hundreds of one-dimensional numerical simulations were conducted by changing the $\lambda$ and $\delta$ of the initial stratification in order to obtain the $\lambda$-$\delta$ regime diagram. The numerical simulation results enabled evaluation of the theory as well as the analysis of the transient behavior of auto-ignition/flame front propagation. The aforementioned framework was applied to different averaged temperatures (830–1000 K) and pressures (20 and 50 bar) for relatively sensitive (PRF-E) and nonsensitive (PRF87) gasoline surrogates. Finally, a diagnostic approach was presented for creating the predictive combustion mode maps for various gasoline surrogates under thermal stratification. Even within the limited scope of the present numerical/theoretical study, the studied portion of the $\lambda$-$\delta$ phase space is noted to pose a rich topological structure. We report findings on two combustion modes which match the previous picture presented in the literature. The following conclusions are made:

1. The theoretical $\beta \approx 1$ curves are consistent with the resulted phase border from 1D simulations for both surrogates. Furthermore, both the theory and simulations depict that higher stratification levels (larger $\lambda$ and $\delta$) increase the prevalence of the deflagration regime compared to the spontaneous ignition mode.

2. Negative temperature coefficient (NTC) in the fuel chemistry increases the dominance of spontaneous ignition regime by reducing the variation of reactivity versus temperature and shifting the phase borders to higher $\lambda$ and $\delta$. Furthermore, as a result of the NTC behavior, the rate of transition between the two combustion modes is decreased, resulting in a broad range of $\lambda$ and $\delta$ that $\beta \approx 1$ applies. The mentioned behavior is also obtained using the theoretical analysis by showing high sensitivity of the $\beta$-curve borderline to the values of $0.4<C_{\beta }<3.0$ coefficient.

3. The area of spontaneous ignition zone in the $\lambda$-$\delta$ regime diagram is used as a qualitative representation of the knock propensity of the fuels. Using the present approach, observed relations between the fuel sensitivity and the knock resistance at different pressures and temperatures are found to be inline with previous studies. Moreover, the predicted knock propensities of the studied fuel mixtures conform with their similar MON and different RON numbers.

4. Estimated front displacement speeds from transient analysis of flame depicted that the closer the conditions are to the spontaneous ignition region, the higher the front propagation speeds are compared to the laminar flame speed, such a situation represents a blended combustion mode where the flame propagation is assisted with ignition (HTC) of reactants upstream of the flame front, indicating an spontaneous ignition assisted deflagration.

5. For different points on the $\lambda$-$\delta$ regime diagram, the transport budget terms analysis illustrates a reaction term dominated (spontaneous ignition) and a balanced reaction-diffusion term dominated propagation mode (deflagration). Further investigations reveals the possibility of a secondary ignition after the initial ignition kernel for both surrogates. This implies a transition from the deflagration mode to the spontaneous ignition mode, through which the unburnt charge ahead of the flame front is consumed by the secondary ignition due to the heating compression of the approaching fronts and high reactivity of the mixture.

Acknowledgments

We would like to acknowledge Aalto University and CSC (Finnish IT Center for Science) for providing the computational resources.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Additional information

Funding

The present study has been financially supported by Neste Corporation, Finland, and the Academy of Finland 18024,332784 (grant numbers 318024 and 332784).

Appendix A

Chemical mechanism validation

In this section, we present the chemical mechanism validation. Here, the utilized chemical mechanisms are chosen based on the comparison of the calculated IDT and the 1D laminar flamea speed with available experimental data (Dirrenberger et al. 2014; Fikri et al. 2008; Gauthier, Davidson, Hanson 2004; Jerzembeck et al. 2009). As mentioned in section 2.3, constant volume 0D homogeneous reactor and 1D premixed flame modules of Cantera are used for calculation of IDT and flame speed, respectively. In order to maintain the numerical feasibility of the simulations, only reduced mechanisms have been examined in the present study. Another important aspect considered in our selection was the validity of the mechanisms at high pressures. Details on the studied mechanisms are provided in the following, first for the PRF (PRF87) and then for the PRF-E mixture.

Regarding to the PRF, Jerzembeck et al. (2009) developed two reduced mechanisms for the PRF combustion, one for only high temperature chemistry (HTC) with 99 species and 669 reactions, and the other for combined HTC-LTC with 203 species and 1071 reactions, based on the Lawrence Livermore (LLNL) detailed mechanisms (Curran et al. 1998, 2002). Moreover, they provided an experimental data set for laminar flame speeds of $n$-heptane, iso-octane and PRF 87 mixtures at engine relevant conditions. Wang, Yao, and Reitz (2013) proposed a reduced mechanism with 73 species and 296 reactions for diesel/gasoline combustion based on the LLNL detailed mechanism (Mehl et al. 2011b). A 171 species and 861 reactions PRF skeletal mechanism was developed by Luong et al. (2013) using the LLNL detailed mechanisms (Curran et al. 1998, 2002) and further validated in Perfectly Stirred Reactors for IDT and extinction residence times. Moreover, Mehl et al. (2011a) provided a reduced version of the detailed LLNL mechanism for CFD engine modeling. The mechanism has been reduced and validated based on a target non-oxygenated gasoline surrogate containing four components.

Figures A1a and Figure A1b show the comparison between the PRF mixture, the calculated IDT values, and measurements (Gauthier, Davidson, Hanson 2004) at 20 and 50 bar, respectively. In addition, to compare the performance of the reduced mechanisms, results from a recent gasoline surrogate detailed mechanism by Fang et al. (2022), comprised of 1543 species and 8725 reactions, are included. According to Figure A1 the examined mechanisms are in a good agreement with the experimental data. However, the NTC behavior has been predicted differently. In fact, Wang mechanism did not predict any NTC behavior for both of the investigated pressures. Note that results of the Jerzembeck HTC mechanism are not provided since the mechanism did not include LTC related reactions.

Figure A1. PRF 87, ignition delay time versus temperature, circles display experimental data. (Gauthier, Davidson, Hanson 2004).

Furthermore, the mechanisms have been evaluated against experimental high pressure flame speed data (Jerzembeck et al. 2009). According to the literature, the steady state solution of premixed flame using solvers such as Cantera would not converge under auto-ignitive conditions (high pressures and intermediate temperatures). The problem is relevant to the low temperature reactions initiating a secondary flame upstream of the primary reaction zone and violating the cold boundary condition assumption for steady state solutions. This problem can be addressed by excluding the LTC reactions from the mechanism for prediction of flame propagation only, e.g. Jerzembeck HTC mechanism. However, recent studies of Sankaran (2015) and Krisman, Hawkes, and Chen (2018) shows that for a sufficiently small domain, the steady state solution is converging to reaction rates associated with single stage ignition. Therefore, in the present study the same approach has been taken to find the correct domain length for convergence.

Figure A2 depicts the predicted laminar flame speeds versus equivalence ratio against experimental data at 20 bar and 298 K initial temperature. It is observed from Figure A2 that the Wang, Luong, and LLNL reduced mechanisms underpredict the laminar flame speed. As stated earlier, except for the Polimi mechanisms all other mechanisms are developed based on the LLNL detailed mechanism, however, only Jerzembeck et al. (2009) considered the high pressure laminar flame speeds at the development stage, and hence their mechanisms results deviate from others.

Figure A2. Laminar flame speed versus equivalence ratio for PRF 87 combustion at 298 K and 20 bar, circles display experimental data. (Jerzembeck et al. 2009).

Figure A2 shows that although Polimi Biogasoline mechanism underpredicts the laminar flame speed for rich mixtures, results of both mechanisms are in agreement with experiments at stoichometric equivalence ratio (considered in the present study). Based on the above validation results, Polimi reduced mechanism has been chosen for the PRF mixture related simulations.

Regarding to the PRF-E, Ranzi et al. (2014), developed a number of skeletal mechanisms for typical surrogate mixtures using reaction flux analysis. They also extended the reaction mechanisms for biogasolines, providing a 171 species and 3751 reactions mechanism including ethanol (POLIMI BIOGASOLINE 171). The mechanism has been validated against alcohols atmospheric laminar flame speeds and the PRF-ethanol ignition delay time measurements (Fikri et al. 2008). Based on the POLIMI 1412 mechanism (Ranzi et al. 2012), Stagni et al. (2016) proposed a reduced mechanism with 115 species and 856 reactions using a novel automatic mechanism reduction method.

Further validations have been conducted in order to assess the Polimi Biogasoline mechanism performance using the PRF-E mixture IDT and pure ethanol laminar flame speed data sets (Dirrenberger et al. 2014; Fikri et al. 2008). Figure A3 shows the comparison between the PRF-E ignition delay data and the experiments. Moreover, the predicted laminar flame speed of pure ethanol against experimental data is shown in Figure A4.

Figure A3. Ignition of PRF-E/air stoichiometric mixtures predicted by Polimi Biogasoline mechanism at different pressures, symbols depict experiments (Fikri et al. 2008).

Figure A4. Laminar flame speeds of pure ethanol predicted by Polimi Biogasoline mechanism, symbols depict experiments (Dirrenberger et al. 2014).

Figure A3 and A4 show that the PRF-E IDTs at different pressures using the Polimi Biogasoline mechanism are in good agreement with experimental data. Based on present validation results, the Polimi Biogasoline mechanism has been chosen for the PRF-E mixture simulations.