Temperature Stratification Induced Ignition Regimes for Gasoline Surrogates at Engine-Relevant Conditions

ABSTRACT

End-gas auto-ignition leading to knocking combustion is one of the major barriers to achieving higher thermal efficiencies in downsized boosted spark-ignition engines. Despite the available framework addressing hotspot-induced ignition (detonation peninsula), a quantitative investigation on hotspot-induced auto-ignition of gasoline surrogates is yet to be done. In particular, the effect of negative temperature coefficient (NTC) chemistry on the distribution of the ignition modes in the detonation peninsula is still missing. Using the established one-dimensional (1D) theoretical and computational framework, the effect of average temperature (including NTC range), initial pressure, and ethanol addition are investigated. Moreover, appearance of NTC chemistry-related events i.e. coolspots, secondary ignition kernels, and off-centered ignition are analyzed using 1D simulations. The results are as follows. 1) NTC chemistry affects the distribution of ignition regimes in detonation peninsula and the dynamics of the front propagation via altering the reactivity gradient. 2) NTC chemistry increases the temperature gradient range associated with the detonation regime. 3) NTC may inhibit detonation development by simultaneously promoting the spontaneous/supersonic ignition modes. 4) An ethanol blend decreases the knock propensity; however, lower ignitability may promote detonation development and the appearance of strong shock waves. 5) Finally, detonation may result in a normal knock at lower initial pressures (20 bar). However, at elevated initial pressures (50 bar), detonation is noted to yield pressure intensities resembling super-knock.

KEYWORDS:

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

Introduction

Recent progress in engine efficiency, as well as CO₂ emission control, has derived the development of Spark Ignition (SI) engines toward downsized and boosted operating scenarios. However, the primary obstacle for achieving high compression ratios is the abnormal combustion phenomenon i.e. the risk of knock and super-knock in normal and turbo-charged engines (Kalghatgi and Bradley 2012; Kalghatgi et al. 2009). Engine knock is associated with a rapid heat release and the subsequent pressure oscillations initiated by auto-ignition of the fuel/air mixture (Amann, Alger, and Mehta 2011). Moreover, intermittent pressure waves and ignition generated by developing detonation or deflagration to detonation transition (DDT) events can lead to super-knock (Wang et al. 2015; Wei et al. 2021). Therefore, in the SI engines context, understanding the characteristics of auto-ignition and its various propagation mechanisms is of utmost importance. In this study various auto-ignition propagation regimes of gasoline surrogates are investigated using a one-dimensional (1D) computational framework under normal and boosted SI engine relevant conditions.

Auto-ignition characteristics (knock propensity) of gasoline blends are measured by Research Octane Number (RON) (Astm 2021) and Motored Octane Number (MON) (Astm 2022) tests. Primary Reference Fuels (PRFs), i.e. binary mixtures of $\mathit{n}$-heptane and iso-octane, are the gasoline surrogates used as a reference (RON = MON) to measure octane qualities of the fuels in these tests. Moreover, ignition characteristics of sensitive (S = RON – MON $\ne$ 0) fuels can be deduced from the octane index (OI) measured under modern engine operational conditions (Bradley and Head 2006; Bradley, Morley, and Walmsley 2004). With relevance to the present work, bio-derived oxygenates (e.g. ethanol) are reported to improve the anti-knock quality of gasoline (Cheng et al. 2020; Pan et al. 2021). Emissions reduction, low energy density, and high octane sensitivity make ethanol a beneficial additive for gasoline (Cheng et al. 2020). In this respect, PRF87 (13$\mathrm{\%}$ $\mathit{n}$-heptane/87$\mathrm{\%}$ iso-octane volume, PRF) and PRF-ethanol (18$\mathrm{\%}$ $\mathit{n}$-heptane/62$\mathrm{\%}$ iso-octane/20$\mathrm{\%}$ ethanol volume, PRF-E) mixtures are chosen as gasoline surrogates in the present work.

In SI engines, auto-ignition can occur in random locations and time instances in the unburnt mixture between the flame and cylinder walls, i.e. the end-gas. The potential inhomogeneity of temperature and mixture composition affects the reactivity distribution in the end-gas and forms localized, isolated exothermic regions or ”hotspots.” Originally outlined by Zeldovich (1980), five modes of propagation were identified by Gu, Emerson, and Bradley (2003) for H₂-CO/air hotspot-induced ignition fronts using 1D simulations. Thermal explosion (spontaneous ignition), supersonic ignition, developing detonation and subsonic ignition are the identified regimes initiated from a hotspot. Using a similar approach, Bradley et al. (2002) defined a quantitative regime diagram (detonation peninsula) based on two dimensionless parameters for H₂-CO/air mixtures under Controlled Auto-Ignition (CAI) engines relevant conditions.

The original detonation peninsula (Bradley et al. 2002) and the associated 1D framework have been used to assess the knock and super-knock tendencies under SI engine relevant conditions (Bates et al. 2016; Bradley and Head 2006; Kalghatgi and Bradley 2012). More recently, Gorbatenko, Bradley, and Tomlin (2021) reported the anti-knock properties of $\mathit{n}$-butanol/air and its blends with Toluene Reference Fuel (TRF) using the original peninsula. Their results indicated that blending with $\mathit{n}$-butanol improves the knock resistance at 50 bar. Moreover, Su, Dai, and Chen (2021) reproduced the detonation peninsula for methane/air mixtures. The detonation development regimes, obtained by three detailed kinetic models, showed qualitative agreements with each other and with the original peninsula. However, they found quantitative discrepancies between the different chemical mechanisms results and the original peninsula. The role of chemical reactivity and energy density on the hotspot-induced reaction wave propagation in stoichiometric mixture of methane/air mixture was studied by Pan et al. (2021). Their results revealed that while high energy density can promote the detonation development, the chemical reactivity has the most significant effect on the induced ignition regimes. The combined effect of energy deposit, acoustic and ignition time scales on the combustion regimes inside the hotpot was studied by Kiverin et al. (2013). Utilizing detailed chemistry and transport models for mixture of H₂/O₂, they demonstrated that the interaction between these parameters determines the gas-dynamic processes, formation and steepness of the temperature gradient, and ignition front speed.

As noted earlier, despite the successful evaluation of the ignition regimes using the 1D framework, it has been shown that the boundaries of the different regimes in the original detonation peninsula are sensitive to the choice of fuel (Dai et al. 2015; Gorbatenko, Bradley, and Tomlin 2021; Pan et al. 2019; Robert et al. 2019; Rudloff et al. 2013) and kinetic mechanism (Liberman et al. 2019; Su, Dai, and Chen 2021). Therefore, it is necessary to reproduce the detonation peninsula to explore the ignition regimes for different fuel-oxidizer mixtures using the 1D framework.

In particular, for large hydrocarbon mixtures with low-temperature chemistry (LTC), multistage ignition, i.e. low-, intermediate- and high-temperature ignition (LTI, ITI and HTI), may complicate reaction-pressure wave interactions (Dai et al. 2015; Ju et al. 2011). Bates et al. (2016) illustrated the detonation peninsula for iso-octane/air and $\mathit{n}$-heptane/air mixtures, however, the results were not fully compared with the original peninsula. Robert et al. (2019) studied different ignition regimes for a TRF/air mixture by using tabulated ignition data at high temperatures. The results indicated that the generated peninsula for the TRF lies within the original peninsula for restricted hotspot radii with values relevant to SI engines. Furthermore, auto-ignition characteristics of $\mathit{n}$-heptane were studied at initial temperatures within and below Negative Temperature Coefficient (NTC) range by Dai and coworkers (Dai and Chen 2015; Dai et al. 2015). In the presence of NTC chemistry, coolspots were identified as sources of ignition fronts leading to knock in engines. Moreover, multiple-ignition kernels/fronts, shock waves, and detonation waves were identified and found to be generated by multistage heat release (MSHR). Nonetheless, Dai and coworkers did not compare the identified ignition regimes for the $\mathit{n}$-heptane/air mixture against the original peninsula.

With relevance to the studied fuel mixtures in the present work, appearance of different chemical length scales originated by $\mathit{n}$-heptane/oxygen/diluent mixture’s MSHR was shown through the Zeldovich-von Neumann-Döering (ZND) detonation simulations by Liang, M´evel, and Law (2018). Recent 1D transient simulations by Zhao, Ren, and Zhang (2021) demonstrated low frequency periodic amplification of the shock waves by onset of DDT in detonation of off-stoichiometric $\mathit{n}$-heptane/air mixtures. Moreover, locally intensified pressure waves were identified close to the walls of a 1D reactor as a result of NTC chemistry for the $\mathit{n}$-heptane/air mixture (Terashima and Koshi 2015). Furthermore, 1D simulations of the $\mathit{n}$-heptane/air mixture (Terashima, Matsugi, and Koshi 2017) revealed that a hotspot with high reactivity can enhance the reactivity gradient in the end-gas leading to strong shockwave generation, while low reactivity hotspot suppresses the pressure wave development, by lowering reactivity gradient in the end-gas. Recently, zero-dimensional (0D) and 1D study of Yu, Zhang, and Dai (2021) showed that water vapor dilution may decrease the NTC effect in hotspot-induced ignition regimes of $\mathit{n}$-heptane/air mixtures. Moreover, a positive correlation was found between the hotspot radius and maximum pressure in detonative cases, indicating strengthened coupling between the reaction and pressure waves. The prevalence of deflagrative versus auto-ignitive combustion regimes of thermally stratified gasoline surrogates (i.e. PRF and PRF-E) including the diffusion effect prior to ignition was studied by Shahanaghi et al. (2022). The proposed theoretical and 1D numerical framework has been used to estimate and compare the knock propensity of studied mixtures. The results indicated that NTC chemistry increases the dominance of the auto-ignition regime, while it promotes a blended auto-ignition assisted deflagration mode.

According to the present literature review, the relative importance of different ignition regimes in the detonation peninsula are sensitive to the thermo-chemical composition of the fuel-air mixture. However, only limited studies reproduced the detonation peninsula for hotspot initiated ignition of large hydrocarbons with the presence of low-temperature chemistry (LTC). In addition, a comparative study on the ignition regime diagrams of practical gasoline surrogates under SI engine relevant conditions is yet to be conducted. Finally, due to the complex combustion chemistry, the knocking mechanisms of gasoline and gasoline-ethanol blends, are not fully understood, and they require further investigations.

Therefore, in order to address the research gap outlined above and in continuation of our previous study (Shahanaghi et al. 2022), here, we investigate hotspot-initiated ignition of gasoline surrogates, under normal and boosted SI engine relevant conditions. Using the established 1D theoretical and numerical framework (Bradley et al. 2002), different pressures, temperatures and gasoline surrogates are examined in a systematic way. Furthermore, the results are analyzed with the help of 0D homogeneous reactor data and 1D ZND data on the detonation structure. The main objectives of the current work are listed as follows. 1) Extend our previous work (Shahanaghi et al. 2022) on the prevalence of deflagration versus auto-ignition in the SI engines context to study various hotspot-induced auto-ignition regimes. 2) Identify the detonation regime for gasoline surrogates (PRF and PRF-ethanol) by using different approaches based on the hotspot pressure. 3) Construct the ignition regime diagrams for the gasoline surrogates under normal and boosted engines relevant conditions. 4) Investigate the effect of NTC chemistry on the prevalence of various ignition modes as well as the appearance of different propagation patterns inside the hotspot. 5) Find the knock propensity and intensity, as well as the responsible mechanisms in the studied cases by comparing the estimated pressure intensity levels and available knock metrics for SI engines.

Methodology

The present study continues the authors previous investigation (Shahanaghi et al. 2022) on the prevalence of deflagration versus auto-ignition combustion modes of thermally stratified gasoline surrogates under SI engine relevant condition. Here, we study the prevalence of auto-ignition propagation modes, i.e. subsonic ignition, detonation and supersonic/spontaneous ignition propagation (Zeldovich 1980) of the gasoline surrogates within a hotspot under the SI engine end-gas related pressures and temperatures. The present work utilizes the theoretical and numerical framework in Refs. (Gu, Emerson, and Bradley 2003; Zeldovich 1980) to identify effective parameters in the hotspot ignition process where hundreds of 1D numerical simulation are conducted to study different hotspot ignition regimes for gasoline surrogates. More discussion is provided for these steps in the following sections.

Theoretical background

According to Zeldovich (1980), ignition front propagation speed, $S_{\mathit{ign}}$, can be deduced from the inverse of ignition delay time (IDT) gradient of the premixed reactants. Assuming a compositionally homogeneous mixture, IDT can be considered as a function of temperature and thereby, $S_{\mathit{ign}}$ can be written as

(1) $S_{\mathit{ign}}=\frac{1}{\left|\mathrm{\nabla }{\tau }_i\right|}={\left(\frac{d{\tau }_i}{\mathit{dT}}\frac{\mathit{dT}}{\mathit{dx}}\right)}^{-1}=(\mathit{\alpha }\nabla \mathit{T}{)}^{-1},$(1)

where ${\tau }_i$ is IDT, $\mathit{\alpha }=\frac{d{\tau }_i}{\mathit{dT}}$ represents the IDT sensitivity to temperature variations and $\mathrm{\nabla T}$ is temperature gradient. Two types of waves are generated at the instance of auto-ignition; a reaction wave representing the propagation of chemical reactions and a pressure wave (shock wave) generated by thermal explosion of the mixture (reaction shocks). According to the Zeldovich theory (Zeldovich 1980), ignition fronts with respect to their propagation speeds can be divided into four groups. First, deflagration which is associated with steep temperature gradients and propagating with the laminar flame speed. Second, subsonic ignition in which the reaction front is faster than deflagration but slower than the pressure wave. Third, supersonic ignition in which the reaction wave propagates faster than the pressure wave. Zeldovich characterized thermal explosion (spontaneous ignition) as a limiting case for supersonic ignition in which the reaction front approaches an infinite propagation speed. Finally, detonation, in which $S_{\mathit{ign}}$ is close to the local speed of sound ($U_a$). In such cases, an inherent coupling between the reaction waves and shock waves may initiate detonation development (Bates et al. 2016; Gu, Emerson, and Bradley 2003). Considering Equation 1(1) $S_{\mathit{ign}}=\frac{1}{\left|\mathrm{\nabla }{\tau }_i\right|}={\left(\frac{d{\tau }_i}{\mathit{dT}}\frac{\mathit{dT}}{\mathit{dx}}\right)}^{-1}=(\mathit{\alpha }\nabla \mathit{T}{)}^{-1},$(1) , the critical temperature gradient ($\mathrm{\nabla }{T}_c$) associated with such a scenario is written as

(2) $\mathrm{\nabla }{T}_c=(\mathit{\alpha }{U}_a{)}^{-1}.$(2)

The original theory by Zeldovich (1980) has been further developed by Bradley and coworkers (Bradley et al. 2002; Gu, Emerson, and Bradley 2003) to study the ignition regimes inside a spherical hotspot with radius $r_0$. A 1D numerical framework as well as a regime diagram known as detonation peninsula was defined based on two dimensionless parameters ($\mathit{\xi }$, $\mathit{\epsilon }$). The original detonation peninsula was developed for H₂-CO/air combustion inside controlled auto-ignition (CAI) engines (Bradley et al. 2002; Gu, Emerson, and Bradley 2003) and later examined for various engine types (Bates et al. 2016). The first dimensionless parameter, $\mathit{\xi }$, indicates the degree of coupling between the acoustic wave and the reaction front, written as

(3) $\mathit{\xi }=\frac{U_a}{S_{\mathit{ign}}}=\frac{\mathrm{\nabla T}}{\mathrm{\nabla }{T}_c}.$(3)

As mentioned earlier, theoretically, detonation kicks in as soon as $U_a=S_{\mathit{ign}}$ or $\mathit{\xi }$ = 1.0. However, due to thermal and species diffusion during the induction period, appearance of detonation is not limited to $\mathit{\xi }$ = 1.0 criterion. In practice, depending on the mixture’s reactivity, a wider range of $\mathit{\xi }$ values can lead to a developing detonation inside a hotspot. Therefore an upper limit (${\xi }_u$) and a lower limit (${\xi }_l$) are defined and the confined zone between these limits indicates the detonation peninsula. Moreover, ${\xi }_u$ (${\xi }_l$) indicates the border between the detonation regime and the subsonic (supersonic) ignition mode.

The second dimensionless parameter, $\mathit{\epsilon }$ indicates the ratio between the chemical timescale i.e. excitation time (${\tau }_e$) and the acoustic time scale $r_0/U_a$ in the hotspot,

(4) $\mathit{ϵ}=\frac{r_0/U_a}{{\tau }_e}.$(4)

Equation 4(4) $\mathit{ϵ}=\frac{r_0/U_a}{{\tau }_e}.$(4) , demonstrates the rate at which chemical energy is transferred to the pressure wave. In practice, a minimum distance (here $r_0$) is required from the ignition onset to the appearance of detonation (detonation run-up distance).

In addition to the mixture’s reactivity (i.e. ${\tau }_e$ and ${\tau }_i$) and its variation versus temperature, recently Pan et al. (2021); Pan et al. (2021) studied the role of fuel energy density ($E_c$) on ignition and knocking initiated by hotspots. Energy density was found to promote the detonation development, and a positive correlation between the detonation peak pressure and energy density was observed. Energy density can be calculated using the ideal-gas assumption and lower heating value of fuels,

(5) $E_c=\frac{n{H}_m}{\mathit{V}}=\frac{{\mathit{P}}_{\mathit{fuel}}{H}_m}{\mathit{T}{R}_u},$(5)

where n is the number of moles, $H_m$ is the lower heating value and $\mathit{V}$ is the volume, $P_{\mathit{fuel}}$ is the partial fuel pressure and $R_u$ = 8.314 J/mol.K is the universal gas constant.

It should be noted that in the above equations, values of ${\tau }_i$, ${\tau }_e$, $U_a$, $\mathit{\alpha }$ and $H_m$ as well as the derived parameters, i.e. $\mathit{\xi }$, $\mathit{\epsilon }$, $S_{\mathit{ign}}$ and $E_c$ correspond to a reference temperature ($T_{\mathit{ref}}$). It is shown that $\mathit{\xi }$ and $\mathit{\epsilon }$ values are sensitive to the choice of $T_{\mathit{ref}}$ (Robert et al. 2019). Therefore, to make the results comparable with the original peninsula presented by Bradley et al. (2002), in the present work, $T_{\mathit{ave}}=\mathit{T}\left(\mathit{r}=r_0/2\right)$     is chosen as $T_{\mathit{ref}}$. In the present study, the ${\tau }_i$ values are calculated using the 0D constant volume reactor model of Cantera (Goodwin et al. 2018). Moreover, the generated ${\tau }_i$ vs temperature data are used to estimate $\mathit{\alpha }$ values. Furthermore, ${\tau }_e$ is considered as the time until reaching 5% of the maximum heat release in the constant volume reactor (Bradley et al. 2002). Finally, $U_a$, and $H_m$ are calculated by using the thermodynamic library of Cantera.

Transient 1D simulations

In this part, first, the initialization and expected propagation patterns in 1D constant volume domain is explained. Then, the solver for 1D simulations is introduced and the solver and chemical mechanism validation details are provided.

Initialization and propagation patterns

The setup in the present paper concerns formation of auto-ignition waves by a hotspot inside a 1D constant volume domain, Figure 1a, representing a spherical combustion chamber. The hotspot is modeled as a localized linear temperature gradient with a specified average temperature ($T_{\mathit{ave}}$), temperature gradient ($\mathrm{\nabla T}$) and radius ($r_0$). The initial pressure and mixture composition are considered to be homogeneous along the 1D domain.

Figure 1. A schematic presentation of the combustion process shows, a) hotspot/coolspot initial T profile, b) centered ignition and regular propagation inside the hotspot/coolspot, c) off-centered primary ignition due to NTC and double ignition front propagation, and d) occurrence of primary ignition and subsequently, formation of the secondary ignition kernel due to NTC.

For mixture composition we select two different gasoline surrogates, same as those in our earlier study (Shahanaghi et al. 2022), in order to address the effect of fuel octane sensitivity (S) on the detonation and knock characteristics. A binary blend of iso-octane/$\mathit{n}$-heptane (87%13% volume, PRF) and a ternary blend of iso-octane/$\mathit{n}$-heptane/ethanol (18%/62%/20% volume, PRF-E) are considered for numerical simulations. It is worth noting that addition of ethanol reduces the NTC behavior in the PRF mixture. This resembles the increase in knock resistance of a commercial gasoline by using oxygenated additives. The selected PRF mixture has a similar MON (87) to that of gasoline while the PRF-E mixture mimics the oxygenated gasoline with MON = 89 and RON = 95.1. It is worth noting that modern engines operate at higher pressures and lower temperatures than MON conditions, which is more similar to RON-like conditions. Therefore, two initial pressures, i.e. 20 and 50 bar are considered for the selected surrogates as the representation of normal and boosted engine scenarios.

The anticipated ignition patterns inside the hotspot after primary ignition are shown in Figure 1b–d. Inside the hotspot, the mixture is initialized with higher reactivity in the hotspot center ($\mathit{r}=0$), Figure 1a; therefore, the initial ignition is expected to occur in the center leading to a reaction front propagating outwards in the radial direction, Figure 1b. However, in practice, such a regular propagation pattern is prone to changes due to complex ignition and burning properties of the studied mixtures. Large hydrocarbons utilized in engines and their surrogates, e.g. PRF, usually exhibit low and intermediate temperature chemistry. Such a characteristic feature can manifest multistage ignition, i.e low and intermediate ignition, as well as NTC behavior.

In the presence of strong NTC chemistry, regions at relatively lower temperature will ignite earlier than the rest of the mixture. In such scenarios and in order to maintain an outwardly propagating pattern, the 1D problem is initialized with a locally positive temperature gradient named as a coolspot, Figure 1a. Nevertheless, since $T_{\mathit{ave}}$, and temperature variation ($\mathrm{\Delta T}$ $=\mathrm{\nabla T}$ $\times$ $r_0$) are pre-determined variables in the present work and considering the presence of NTC, the localized high reactivity region may shift to an arbitrary location alongside the radius. In such scenarios, an”off-centered ignition” kernel may appear inside the hot/cool spot leading to the formation of right and left propagating reaction waves, as illustrated in Figure 1c. In addition, under relatively high-temperature stratification levels, the NTC chemistry may result in non-monotonic reactivity stratification inside the hotspot. In such occasions, the primary ignition occurs at the center of hotspot. However, later, due to the accumulation of low and intermediate heat release, a secondary ignition kernel may appear inside the hotspot leading to a complex system of reaction fronts propagating inwards or outwards, Figure 1d. More details on the selected initial pressures and average temperatures for 1D simulations using the 0D analysis are provided in Section 3.1.

Numerical methods and governing equations

In the present paper, an open source computational fluid dynamics (CFD) package OpenFOAM (Weller et al. 1998) is used to solve the 1D hotspot simulation problem. The solver in the present paper is a fully compressible, density-based flow solver with explicit time stepping based on the rhoCentralFoam (Greenshields et al. 2010) solver. As a part of OpenFOAM, rhoCentralFoam is a density-based solver utilizing second-order central, Kurganov and Tadmor (Kurganov and Tadmor 2000), discretization scheme for capturing shock waves in high speed non-reacting flows. Here, the time integration of rhoCentralFoam is modified to the second-order explicit Runge Kutta scheme. Moreover, similar to Ref. (Khodadadi Azadboni, Heidari, and Wen 2018), the species transport equations together with chemical reaction source terms are included in the solver to be able to solve reacting detonation problems. The chemical source terms are calculated using direct integration of the finite rate chemistry. The Jacobian matrix for the system of kinetic ODEs is calculated analytically using pyJac (Niemeyer, Curtis, and Sung 2017) and its recent implementation to the OpenFOAM framework (Morev et al. 2022). These numerical techniques have been successfully used for various compressible reactive flow simulations (Khodadadi Azadboni, Heidari, and Wen 2018; Morev et al. 2022; Zhao, Ren, and Zhang 2021). The formulation of the flow field governing equations, i.e. continuity, momentum, species concentration and enthalpy are as follows:

(6) $\frac{\mathrm{\partial \rho }}{\mathrm{\partial t}}+\mathrm{\nabla }\cdot (\mathit{\rho u})=0,$(6)

(7) $\frac{\mathrm{\partial }(\mathit{\rho u})}{\mathrm{\partial t}}+\mathrm{\nabla }\cdot (\mathit{\rho uu})=\mathrm{\nabla }\cdot \mathit{\tau }-\mathrm{\nabla p},$(7)

(8) $\frac{\mathrm{\partial }(\rho Y_k)}{\mathrm{\partial t}}+\mathrm{\nabla }\cdot (\mathit{\rho u}{Y}_k)=\mathrm{\nabla }\cdot (\mathit{\rho D}\nabla Y_k)+{\dot{\omega }}_{\mathit{k}},$(8)

(9) $\frac{\mathrm{\partial }(\mathit{\rho h})}{\mathrm{\partial t}}+\mathrm{\nabla }\cdot (\mathit{\rho uh})-\frac{\mathrm{\partial p}}{\mathrm{\partial t}}=\mathrm{\nabla }\cdot (\mathit{\rho }{\alpha }_t\mathrm{\nabla h})+{\dot{\omega }}_{\mathit{h}}.$(9)

In the above equations, $\mathit{\rho }$, $\mathit{u}$, $\mathit{p}$, $Y_k$, $\mathit{h}$, $\mathit{\tau }$, ${\alpha }_t$, and $\mathit{D}$ denote density, velocity, pressure, mass fraction of the species $\mathrm{k}$, sensible enthalpy, viscous stress tensor, mass and thermal diffusivities, respectively. The production rates of each species and heat release rate (HRR) are represented by ${\dot{\omega }}_{\mathit{k}}$ and ${\dot{\omega }}_{\mathit{h}}$, where ${\dot{\omega }}_{\mathit{h}}$=${\mathrm{\Sigma }}_k\mathrm{\Delta }{h}_{f,k}^0{\dot{\omega }}_{\mathit{k}}$, in which $\mathrm{\Delta }{h}_{f,k}^0$ is the enthalpy of formation. In the present work, estimated Lewis number for fuel/oxidizer mixture is about unity. Therefore, the unity Lewis number assumption ($\mathit{D}={\alpha }_t$) is incorporated for all species. In addition, viscosity and thermal conductivity are calculated with Sutherland law and Eucken approximation (Poling, Prausnitz, and O’Connell 2001). Finally, the ideal gas law and thermal equation of state close the system of equations.

Case setup and validation

The simulation setup for the hotspot-induced ignition transient simulations is a 1D constant volume chamber. A symmetry boundary condition (BC) is applied to the left boundary while wall condition is considered for the right BC. The 1D domain length (L) is considered sufficiently long (L = 4 cm) to avoid appearance of reflected shock waves from the right BC during combustion inside the hotspot (${\mathrm{r}}_0\stackrel{<}{\sim \phantom{{ }_x}}$ 3 cm). The hotspot interior is resolved with uniform grid spacing of (0.75 $\mathit{\mu }$m) while the rest of the domain is discretized with a coarse grid (0.2 mm). The coarse grid resolution outside the hotspot region will damp the small waves reflected from the right BC. Moreover, the utilized grid spacing ensures more than 10 cells inside the Half Reaction Length (HRL) for the ZND structure of detonation. HRL, $L_2$ in this study, is considered as the distance from the leading shock to the location of maximum heat release (Lee 2008) and calculated with SDToolbox (Kao and Shepherd 2008). More discussions regarding the ZND detonation structure of the studied mixtures is provided later in Section 3.2 and Appendix B. Nevertheless, the grid sensitivity studies, provided in Appendix A, showed convergence of the results for the selected grid resolution. In addition, the solver validation for a non-reacting shocktube (Sod 1978) problem as well as a CH₄ hotspot detonation case (Pan et al. 2019) are presented and compared with available data from the literature in Appendix A.

In order to conduct the 1D simulations and by considering the limitation of computational time, two reduced mechanisms are selected for PRF and PRF-E mixtures, respectively. A 115 species skeletal TPRF mechanism presented by Stagni et al. (2016) for PRF and a 171 species skeletal mechanism (Ranzi et al. 2014) for PRF-E are selected. It is notable that in our previous study (Shahanaghi et al. 2022) we validated these mechanisms against experimental data and several kinetic mechanisms including a detailed mechanism by using high-pressure 0D constant volume reactor and 1D laminar flame simulations.

A-Priori analyses

Zero dimensional homogeneous ignition

As noted earlier, in this study, constant volume homogeneous reactor calculation is needed to derive the chemical timescales, i.e. ${\tau }_i$ and ${\tau }_e$, related to the 1D framework. In addition, the chemical composition effect on the auto-ignition characteristics of the gasoline surrogates can be attributed to the constant volume ignition delay time (Sarathy, Farooq, and Kalghatgi 2018). Therefore, 0D simulations are utilized in order to choose the corresponding $T_{\mathit{ave}}$ for 1D simulations under normal and boosted SI engine-related pressures (i.e. P = 20 and 50 bar). The IDT profiles from 0D homogeneous reactor simulations are plotted against temperature for PRF and PRF-E in Figure 2.

Figure 2. 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.

Figure 2a depicts that both surrogates exhibit NTC behavior at lower pressure (20 bar), although PRF shows a higher degree of NTC than PRF-E. However, Figure 2b depicts that at 50 bar, the IDT profile for PRF shows NTC behavior while this behavior is milder for PRE-E. The similarity observed in IDT profiles between the two fuels at low pressure matches with their similar MON number values while different profiles at high pressure indicate the distinction of the mixtures RON values.

Furthermore, five different points (case I-V) are selected in IDT-T diagram for further examination in hotspot detonation simulations. Case I and III ($T_{\mathit{ave}}=1000$ K) are located at relatively high-temperature regions where, the IDT values are similar between the two fuels at lower and higher pressures, respectively. Case II ($T_{\mathit{ave}}=830$ K, P = 20 bar) is located at the NTC region for both PRF and PRF-E mixtures. Moreover, case IV ($T_{\mathit{ave}}=900$ K, P = 50 bar) is located at the NTC region for the PRF mixture while no NTC appears in the PRF-E mixture, Figure 2b. Finally, case V ($T_{\mathit{ave}}=830$ K, P = 50 bar) is located at low-temperature zone for both fuels. Table 1 presents the initial pressure (P), average temperature ($T_{\mathit{ave}}$), IDT (${\tau }_i$), excitation time (${\tau }_e$), IDT gradient versus temperature ($\mathit{\alpha }$), speed of sound ($U_a$) and constant volume pressure ($P_{\mathit{cv}}$) for the selected points from 0D simulations.

Table 1. Conditions for the five cases: pressure, average temperature, IDT, excitation time, IDT gradient versus temperature, acoustic speed, and constant volume pressure.

	Fuel	P[bar]	${\mathbf{T}}_{\mathit{ave}}$[K]	${\tau }_{\mathit{i}}$[ms]	${\tau }_{\mathit{e}}$[$\mathit{\mu }$s]	$\mathit{\alpha }$[s/K]	${\mathbf{U}}_a$[m/s]	${\mathbf{P}}_{\mathit{cv}}$[bar]
Case I	PRF	20	1000	2.9	9	4.6e-5	589.23	67.44
	PRF-E	20	1000	2.2	13	3.5e-5	589.61	67.07
Case II	PRF	20	830	6.7	11	2.2e-5	539.98	78.86
	PRF-E	20	830	11	11	1.4e-5	540.36	78.66
Case III	PRF	50	1000	1.2	7	1.3e-5	589.23	170.71
	PRF-E	50	1000	0.9	8	1.2e-5	589.61	170.30
Case IV	PRF	50	900	1.6	9	4.5e-6	560.75	186.79
	PRF-E	50	900	2.2	10	4.7e-6	561.13	186.11
Case V	PRF	50	830	1.5	8	2.5e-6	539.98	199.56
	PRF-E	50	830	2.3	11	5.0e-6	540.36	199.68

Considering the domain length (L = 4 cm) and acoustic speed, Table 1, acoustic pass-by time scale (${\tau }_a=$ $2\mathit{L}/U_a$) can be estimated $\approx 1.0e^{-4}$ s. From Table 1, the ${\tau }_e\ll {\tau }_a\ll {\tau }_i$ condition indicates that acoustic waves equalize the pressure levels inside the hotspot (2 mm $\stackrel{<}{\sim \phantom{{ }_x}}$ r${ }_0\stackrel{<}{\sim \phantom{{ }_x}}$ 30 mm) prior to ignition, and heat release occurs at approximately constant pressure. However, after the ignition, the heat release occurs at a nearly constant volume process and the temperature increase is followed by a rapid pressure rise. We note that no external energy deposit is considered in the present work. Moreover, given the stratification length scales, the thermal diffusion is not active prior to the ignition (Shahanaghi et al. 2022). Therefore, interactions between the external heat source (Kiverin et al. 2013), thermal diffusion, and ignition are considered to be insignificant in the simulations and ignition front propagates via the Zeldovich gradient theory.

ZND planar detonation profiles

Zeldovich-von Neumann-Doering (ZND) model considers a 1D steady-state detonation propagation. Despite the simplifying assumptions, this model can provide a valuable insight into the chemistry behind the detonation leading shock and its effect on the dynamics of detonation propagation (Zhang 2012). In the present work, the ZND detonation structures are calculated, using the SDToolbox (Kao and Shepherd 2008), for the selected fuels at the specified case conditions in Table 1. ZND heat release rate and temperature profiles of PRF and PRF-E are depicted in Appendix B: Figure b1 versus the distance from the leading shock for the studied cases. Discussions on MSHR, induction length ratio and its impact on dynamics of detonation propagation are provided in Appendix B. Table 2 presents the calculated values of $E_c$, $P_{\mathit{znd}}$, detonation equilibrium pressure ($P_{\mathit{cj}}$), Chapman Jouguet speed ($U_{\mathit{cj}}$), and induction length ratio ($L_1/L_2$) for the studied cases.

Table 2. Conditions for the five cases: pressure, average temperature, energy density, von-Neumann pressure, detonation equilibrium pressure, CJ velocity, and induction length ratio.

	Fuel	P[bar]	${\mathbf{T}}_{\mathit{ave}}$[K]	${\mathbf{E}}_c$[MJ/m^3]	${\mathbf{P}}_{\mathit{znd}}$[bar]	${\mathbf{P}}_{\mathit{cj}}$[bar]	${\mathbf{U}}_{\mathit{cj}}$[m/s]	L1/L2
Case I	PRF	20	1000	19.52	219.5	123.1	1841.74	0.033
	PRF-E	20	1000	19.45	218.9	122.5	1840.97	0.039
Case II	PRF	20	830	23.36	263.6	145.8	1841.70	0.031
	PRF-E	20	830	23.28	262.9	145.5	1840.98	0.041
Case III	PRF	50	1000	48.80	560.3	312.1	1862.07	0.039
	PRF-E	50	1000	48.63	558.7	311.4	1861.20	0.050
Case IV	PRF	50	900	54.00	620.6	343.6	1860.65	0.045
	PRF-E	50	900	53.82	618.9	342.9	1859.81	0.060
Case V	PRF	50	830	58.40	671.6	370.1	1859.80	0.055
	PRF-E	50	830	58.20	669.8	369.3	1858.97	0.077

The E_c values calculated for the PRF-E mixture in Table 2 are lower compared to those of the PRF mixture in case I-V. As a result, a similar trend is evident between the ZND detonation equilibrium ($P_{\mathit{cj}}$) and peak ($P_{\mathit{znd}}$) pressures of the PRF-E and PRF mixtures. Such a trend indicates that the PRF-E mixture is potentially less prone to generating intense pressure waves and, subsequently, knocking in SI engines.

The values of $P_{\mathit{znd}}$ and $P_{\mathit{cj}}$ will be utilized in Section 4.2 to point out the detonation regime in 1D hotspot-induced auto-ignition simulations. Moreover, the presence of unsteady events such as DDT or pulsating detonation will be identified by comparing the hotspot detonation peak pressure and $P_{\mathit{znd}}$. In addition, in Section 4.3, the deviation between $P_{\mathit{znd}}$ and $P_{\mathit{cv}}$ will be used to estimate the pressure intensity in end-gas auto-ignition scenarios leading to detonation.

Results and discussions

Hotspot-Induced ignition regimes

In this section, ignition regime diagrams are provided for 1D modeled hot-spots for case (I-V) in Table 1. We use two methods for the analysis within the hotspot, which are described in the next Section 4.1.1. Figures 3–7 depict the regime diagrams identified using the first method as well as ignition delay profiles for PRF and PRF-E at the specified conditions. As it was discussed earlier in Section 2.1, $\mathit{\epsilon }$ and $\mathit{\xi }$ demonstrate the reactivity and coupling factors in the present simulations. The regime diagrams, Figures 3 and 7(a,b), illustrate the supersonic/spontaneous ignition (indicated by letter P), detonation (red peninsula), and subsonic ignition (indicated by letter B) regimes. Each of the mentioned plots include data from several ($\stackrel{>}{\approx \phantom{{ }_2}}50$) 1D simulations by varying $\mathit{\xi }$ and $\mathit{\epsilon }$. Moreover, the secondary horizontal and vertical axes are conversions of $\mathit{\xi }$ and $\mathit{\epsilon }$ to their corresponding temperature gradient ($\mathrm{\nabla T}$) and hotspot radius ($r_0$) ranges. In the present work, the considered range $0<\mathit{\epsilon }<40$, is chosen to make the regime diagrams consistent and comparable with the original peninsula (Bradley et al. 2002). In addition, a realistic condition at $\mathrm{\nabla T}$ = 2 K/mm and $r_0$ = 5 mm is indicated in Figures 3 and 7(a,b) by the letter N for knocking hotspot ignition in SI engines (Kalghatgi et al. 2009). More details are provided later in Section 4.3.

Figure 3. Ignition regime diagram of case I ($\mathit{P}$ = 20 bar and $T_{\mathit{ave}}$ = 1000 K), red zone indicates the detonation peninsula of a) PRF and b) PRF-E mixtures, respectively. Colored circles display the hotspot pressure calculated for each 1D simulation. Letters P, B, and N indicate supersonic ignition, subusonic ignition and the nominal condition for knock initiation inside SI engines (Kalghatgi et al. 2009). c) IDT distribution at the initial time versus T of PRF and PRF-E mixtures.

Figure 4. Ignition regime diagram of case II ($\mathit{P}$ = 20 bar and $T_{\mathit{ave}}$ = 830 K), red zone indicates the detonation peninsula of a) PRF and b) PRF-E mixtures, respectively. Colored circles display the hotspot pressure calculated for each 1D simulation. Letters P, B, and N indicate supersonic ignition, subusonic ignition and the nominal condition for knock initiation inside SI engines (Kalghatgi et al. 2009). c) IDT distribution at the initial time versus T of PRF and PRF-E mixtures.

Figure 5. Ignition regime diagram of case III ($\mathit{P}$ = 50 bar and $T_{\mathit{ave}}$ = 1000 K), red zone indicates the detonation peninsula of a) PRF and b) PRF-E mixtures, respectively. Colored circles display the hotspot pressure calculated for each 1D simulation. Letters P, B, and N indicate supersonic ignition, subusonic ignition and the nominal condition for knock initiation inside SI engines (Kalghatgi et al. 2009). c) IDT distribution at the initial time versus T of PRF and PRF-E mixtures.

Figure 6. Ignition regime diagram of case I ($\mathit{P}$ = 50 bar and $T_{\mathit{ave}}$ = 900 K), red zone indicates the detonation peninsula of a) PRF and b) PRF-E mixtures, respectively. Colored circles display the hotspot pressure calculated for each 1D simulation. Letters P, B, and N indicate supersonic ignition, subusonic ignition and the nominal condition for knock initiation inside SI engines (Kalghatgi et al. 2009). c) IDT distribution at the initial time versus T of PRF and PRF-E mixtures.

Figure 7. Ignition regime diagram of case V ($\mathit{P}$ = 50 bar and $T_{\mathit{ave}}$ = 830 K), red zone indicates the detonation peninsula of a) PRF and b) PRF-E mixtures, respectively. Colored circles display the hotspot pressure calculated for each 1D simulation. Letters P, B, and N indicate supersonic ignition, subusonic ignition and the nominal condition for knock initiation inside SI engines (Kalghatgi et al. 2009). c) IDT distribution at the initial time versus T of PRF and PRF-E mixtures.

Description of the selected approaches to determine detonation

With regard to the literature, two different identification approaches are used in order to differentiate the detonation regime from the other modes. The identification approaches utilize the information from 0D and ZND simulations to distinguish between the hotspot-induced detonation regime and the other ignition modes in the following sections.

The first method uses weighted average of the maximum pressure inside the hotspot ($P_h$) within the time period between the initial ignition until the combustion completion. Details on calculating the weighted average is provided in Appendix C. A critical temperature value (1800 K) is considered to identify the initial ignition instance in the simulations. Moreover, the combustion completion is defined as when the minimum temperature inside the hotspot reaches the 1D premixed flame products temperature under the same initial conditions ($\approx$ 2600 K). In order to distinguish between the detonative propagation mode and the other combustion modes, the calculated pressures are normalized with the constant volume pressure ($P_{\mathit{cv}}$) values reported in Table 1. The mean value between the maximum and minimum of the normalized pressures, among all 1D simulations in each case, is considered as the threshold to distinguish between the detonation and the other regimes. It is noted that the mean pressure value is $\approx$ 1.5 $P_{\mathit{cv}}$ in all cases. The regions with a normalized pressure value above the mentioned threshold are considered to be inside the detonation peninsula. It is noteworthy that previously, a similar threshold value (i.e. 1.1 $P_{\mathit{cv}}$) was used by Luong et al. (2021) to quantify the heat release associated to detonation. It was shown that the correlation found between the detonation heat release and $\mathit{\xi }$ is insensitive to the choice of the pressure threshold value (Luong et al. 2021).

It is noted that the distance traveled by the ignition front is an effective factor on the final estimated pressure for the hotspots. Therefore $P_h$ represents both time and spatial characteristics of pressure development inside the hotspot. Moreover, the effective distance reflects the occurrence of unsteady scenarios, such as transition of the ignition front to detonation or decaying detonation.

The second method uses the maximum value of pressure ($P_{\mathit{max}}$) during the combustion inside the hotspot as the threshold for identification of the detonation regime. This method is used in Section 4.2. The data is initially filtered and the outliers are detected and removed using interquartile range (IQR) i.e. $\pm 1.5$ IQR (Larose 2005). Furthermore, the $P_{\mathit{max}}$ values are normalized by their respective $P_{\mathit{znd}}$ value in Table 2. Finally, the points on the ($\mathit{\xi }$,$\mathit{\epsilon }$) map wherein the maximum pressure is higher than $P_{\mathit{cj}}$, c.f. Table 2, are considered to be inside the detonation zone (P${ }_{\mathit{max}}/$P${ }_{\mathit{znd}}>$P${ }_{\mathit{cj}}/$P${ }_{\mathit{znd}}\approx 0.55$). A similar criterion has been used in previous studies (Pan et al. 2019, 2021; Su, Dai, and Chen 2021) considering P${ }_{\mathit{max}}>$$P_{\mathit{cj}}$ and the reaction front speeds close to the CJ detonation speed ($U_{\mathit{cj}}$) were considered as the detonation region. It is noted that in the present work, comparison between the ignition front propagation speed and $U_{\mathit{cj}}$ is not utilized as a metric. In particular, the appearance of multiple ignition fronts inside a hotspot due to the NTC chemistry as well as the transition between different ignition regimes are the main reasons that make the calculation of a unique propagation speed inside a hotspot with a specific ($\mathit{\xi }$,$\mathit{\epsilon }$) infeasible. Moreover, previous studies (Pan et al. 2019; Yu and Chen 2015) have noted that in hotspot-induced detonation, the mixture ahead of the leading shock wave undergoes chemical reactions altering the detonation peak pressure and hence, the propagation speed differs from the corresponding $U_{\mathit{cj}}$.

Regime diagram: case i

Here, we look at the obtained detonation regime diagrams and IDT profiles for case I, i.e. T_ave = 1000 K and P = 20 bar. Figures 3(a,b) depict the distribution of different ignition regimes in ($\mathit{\xi }$,$\mathit{\epsilon }$) space. In 1D simulations, the maximum value of $\mathit{\epsilon }$ is limited to 40, resulting to the maximum hotspot size of 2 and 3 cm for PRF and PRF-E, respectively. The difference in the PRF and PRF-E hotspot sizes is due to the longer excitation time of PRF-E in comparison to PRF (Table 1). The dashed lines in Figure 3(a,b) represent the original peninsula from (Bradley et al. 2002) presented for H₂-CO/air mixtures under CAI engine conditions. The red colored peninsula shows the estimated detonation region in the present work using the first identification method. Figure 3(a,b) show that the estimated peninsula lies within the original peninsula’s boundaries for both PRF and PRF-E mixtures. The regions located below and above the detonation peninsula represent the supersonic (P) and subsonic (B) ignition front propagation modes, respectively.

From Figure 3(a,b) it is seen that the considered range of $\mathit{\xi }$ <20 for the simulations includes all the ignition regimes and the size of the detonation zone is similar for the two mixtures. However, Figure 3(a,b) show that $\mathrm{\nabla T}$ related to detonation for PRF (0.1–0.6) is smaller than that for PRF-E (0.2–1) considering the same range of $\mathit{\xi }$. Figure 3c depicts that at T = 1000 K, IDT gradient ($\mathit{\alpha }$) is locally higher for PRF compared to PRF-E mixture. Higher $\mathit{\alpha }$ values result in lower $\mathrm{\nabla }{T}_c$, (Equation 2(2) $\mathrm{\nabla }{T}_c=(\mathit{\alpha }{U}_a{)}^{-1}.$(2) ), and subsequently, lower $\mathrm{\nabla }\mathit{T}$ values for the PRF mixture.

The colored circles in Figure 3(a,b) represent the weighted averaged pressure, i.e. $P_h$, at each simulation point. It is observed that the minimum pressure for both fuels, attributed to the subsonic ignition propagation, is less than $P_{\mathit{cv}}$ ($\approx 67$ bar), c.f. Table 1. Moreover, pressure at the detonation zone may obtain higher values than $P_{\mathit{cj}}$ bar) for both fuels while it is less than P${ }_{\mathit{znd}}(\approx 220$ bar), see Table 2. This is due to the fact that the $P_h$ values are decreased by averaging. The maximum pressure distribution will be discussed in Section 4.2.

Finally, Figure 3(a,b) show that $P_h$ inside the detonation peninsula increases with the hotspot radius, $r_0$. The larger the hotspot size, the higher the chance for coupling between the reaction and the shock fronts, and thereby, the higher probability of detonation. Under a favorable temperature gradient and after a certain $\mathit{\epsilon }$ limit, a detonation front will form inside the hotspot. Therefore, the weighted average pressure $P_h$, increases by increasing the size of the hotspot. Moreover, Figure 3(a,b) show that the maximum pressure is higher for PRF-E compared to PRF despite considering the same hotspot size. Therefore, the higher pressures are due to the larger hotspot sizes in the PRF-E cases compared to the PRF cases at the points with same $\mathit{\epsilon }$.

Regime diagram: case ii

In this section, the regime diagrams for the PRF and PRF-E mixtures at $T_{\mathit{ave}}$ = 830 K and P = 20 bar are studied. Figure 4c shows that $T_{\mathit{ave}}$ is located in the NTC zone for both mixtures. For the PRF mixture, the NTC chemistry results in an off-centered ignition, seen schematically in Figure 1. The cases that show such a behavior are located inside the zone bounded by the dash-dotted line in Figure 4a. Another effect attributed to the NTC region is an increase in the supersonic/spontaneous ignition zone, compared to the original peninsula (Bradley et al. 2002), located at the bottom of the detonation peninsula. This is due to the reduction of $\mathit{\alpha }$ and the shift of ${\xi }_l$ to higher values (note the difference in y-axis scales in Figure 4(a,b)). For the PRF-E mixture, the former effect, i.e. off-centered ignition, did not appear in the ignition profiles while the latter effect is observable in Figure 4b, where the supersonic/spontaneous region is slightly expanded compared to the original peninsula and Figure 3b.

Figures 3 and 4 show that the maximum pressure for both mixtures is higher in case II compared to case I, consistent with a similar trend for $P_{\mathit{znd}}$ and $P_{\mathit{cj}}$ in Table 2. We note that in Figure 4a, the hotspot pressure changes rather smoothly between the detonation zone and the supersonic or subsonic regimes and there is no clear boundary separating these modes.

It is observed that the detonation propensity for the ethanol blend is much lower compared to the PRF mixture. However, for the PRF-E mixture the average pressure reaches higher values inside the detonation peninsula. This can be attributed to the larger hotspot size of the PRF-E mixture (note the difference in secondary x-axis values in Figure 4(a,b)). Moreover, for the PRF mixture, the traveled distance by the detonation wave is reduced by the off-centered ignition phenomenon.

The pressure distribution inside the detonation zone is not uniform as observed from Figure 4a. In particular, inside the detonation region and close to the off-centered ignition border, the pressure decreases. This is due to the fact that NTC increases the ignitability; however, it decreases the effective distance for detonation development as well. It is observed that close to the borderline, a broader area of the mixture is consumed by spontaneous ignition inside a hotspot that is promoted by the NTC chemistry. However, as depicted in Figure 4, increasing the hotspot size and temperature gradient elevates the hotspot pressure. This is due to the appearance of two detonation waves propagating inward and outward of the hotspot. More details on the structure of propagating waves inside the hotspot is provided in Section 4.4 for case IV which can be extended qualitatively for the observed trends in this and other case studies.

Regime diagram: case iii

In this section, the regime diagrams for the PRF and PRF-E mixtures at $T_{\mathit{ave}}$ = 1000 K and P = 50 bar are demonstrated. Figure 5 shows that the estimated detonation peninsula is relatively similar between the two mixtures, similar to case I. However, compared to case I, the pressure inside the detonation peninsula is much higher due to the increased hotspot size and hence, the longer traveling distance of the detonation wave.

Figure 5c shows that IDT versus T for both fuels exhibits a monotonous profile around $T_{\mathit{ave}}$ = 1000. Although the IDT profile of the PRF mixture shows NTC behavior near the T_ave =  900 K, no NTC effect was observed in 1D simulations. The main reasons for the absence of NTC effect in 1D simulations are high $\mathit{\alpha }$ values and the difference in the combustion process of 1D hotspot-induced ignition and 0D homogeneous ignition, wherein transient versus equilibrium processes are involved.

It is observable in Figure 5c that PRF-E has a lower reactivity gradient, $\mathit{\alpha }$, compared to PRF in the range of 950 K < T < 1200 K. Therefore, the detonation peninsula is extended to higher $\mathit{\xi }$ values in PRF-E compared to PRF at 25 $<\mathit{\epsilon }<$ 35, Figure 5(a,b). However, temperature gradients, $\mathrm{\nabla T}$, relevant to the detonation peninsula are slightly larger for PRF compared to PRF-E mixture, Figure 5(a,b). Overall, the estimated peninsula is located within the original peninsula (Bradley et al. 2002) and the ${\xi }_l$ is consistent between the two results.

From Table 2 it is observed that the PRF-E mixture has lower energy density than PRF while PRF-E shows higher reactivity close to $T_{\mathit{ave}}$ as seen from Figure 5. Therefore, higher pressure levels of the PRF-E mixture compared to PRF indicates that in the current case, the positive impact of higher reactivity (lower ${\tau }_i$ values) dominates the negative effects of lower energy density and higher ${\tau }_e$ on the pressure levels.

Regime diagram: case iv

In this section, regime diagrams of the PRF and PRF-E mixtures at $T_{\mathit{ave}}$ = 900 K and P = 50 bar are studied. IDT profiles in Figure 6c show that T_ave is located in the NTC zone for PRF while PRF-E exhibits only a slight inflection. Figure 6(a,b) depict that the detonation peninsula for both mixtures lies within the original curve (Bradley et al. 2002). However, the size of peninsulas differs considerably between the two mixtures.

According to 1D simulations for the PRF mixture, NTC behavior results in the formation of coolspots, where ignition starts from lower temperature and propagates toward high-temperature regions inside the hotspot. The range of $\mathit{\xi }$($\mathrm{\nabla T}$) and $\mathit{\epsilon }$ ($r_0$) attributed to coolspots, indicated by the green dash-dotted line in Figure 6, shows that it includes both supersonic/spontaneous ignition and detonation regimes. The subsonic-ignition regime is located in relatively high values of $\mathrm{\nabla T}$ and $r_0$. However, from Figure 6c, NTC is attributed to a limited temperature range and yet lower $\mathrm{\nabla T}$ values. In addition, off-centered ignition, caused by NTC behavior, appears within the region indicated by the purple dash-dotted line. As depicted in Figure 6a, the off-centered ignition is first appearing in coolspots with low $\mathrm{\nabla T}$ and $r_0$ values, and it is later extended to hotspots and higher values of $\mathrm{\nabla T}$ and $r_0$. Moreover, the NTC effect for the PRF mixture is observed as the appearance of a secondary ignition kernel inside the hotspots. Such a phenomenon will change dynamics of the ignition fronts and later the evolution of pressure inside the hotspot. An example of such scenarios is provided in Section 4.4.5.

Similar to case II, the presence of NTC and the subsequent low-temperature ignition (LTI) boosts the super-sonic/spontaneous ignition regime and hence, ${\xi }_l$ is shifted toward higher values for the PRF mixture. Such an effect is not seen for the PRF-E mixture in Figure 6. In addition, the enhanced reactivity due to NTC extends the upper detonation limit to higher values of $\mathit{\xi }$ for the PRF mixture. Overall, according to Figure 6, the PRF mixture shows a larger detonation zone, 4 $<\mathit{\xi }<$ 35, compared to PRF-E, 1 $<\mathit{\xi }<$ 12.

In the present case, Figure 6a, a large portion of the detonation peninsula is located in the regions where NTC chemistry is active (coolspots, off-centered ignition and secondary ignition). This demonstrates the importance of low and intermediate heat release on the detonation behavior of the PRF mixture under elevated pressure and intermediate temperatures. According to Figure 6a, pressure obtains relatively low values in the detonation regions initiated from coolspots. The low pressure values are due to the reduced effective distance required for detonation development caused by: 1) the increased size of the ignition kernel, due to the decreased $\mathit{\alpha }$ values by NTC chemistry, 2) the occurrence of off-centered ignition and 3) the relatively smaller hotspot radii.

A comparison between the P_h values inside the detonation peninsula for the PRF and the PRF-E mixtures shows higher pressures for the latter while the hotspot radii are similar for the two mixtures. The reduced effective distance for the traveling waves in the PRF mixture due to the off-centered ignition and its correlation with the detonation pressure is the main reason for such a lower pressure.

Regime diagram: case v

Regime diagrams as well as IDT vs T distributions for the PRF and PRF-E mixtures at $T_{\mathit{ave}}$ = 830 K and P = 50 bar are presented in Figure 7. $T_{\mathit{ave}}$ is located in the low-temperature region, Figure 7c, for both mixtures. However, considering case IV ($T_{\mathit{ave}}$ = 900 K), we note that the PRF mixture exhibits NTC behavior in the proximity of $T_{\mathit{ave}}$ in the present case. The IDT profiles vary substantially at temperatures less than $T_{\mathit{ave}}$ while the IDT gradient is small at higher temperatures. Moreover, both IDT values and its gradient are lower for the PRF mixture compared to the PRF-E mixture.

Figure 7(a,b) show that the detonation peninsula is generally located within the boundaries of the original peninsula (Bradley et al. 2002) for both fuels. However, Figure 7a depicts that the predicted peninsula falls below the original line for the PRF mixture (${\xi }_l\approx 0.5$). Such a behavior indicates that under the present initial conditions, the PRF mixture may detonate under lower coupling factors than the critical value (${\xi }_u<1$). Nevertheless, the supersonic/spontaneous ignition regime is limited to a narrow zone compared to the other ignition modes for both fuel mixtures. From Figure 7, considering the ($\mathit{\xi }$,$\mathit{\epsilon }$) coordinates, the upper branch of the detonation peninsula obtains similar values for both fuels (${\xi }_u\approx 25$). Similar to the supersonic mode, the detonation zone is also limited to a narrow region at smaller values of $\mathit{\epsilon }$ while ${\xi }_u$ rapidly increases at $\mathit{\epsilon }\approx 16$ for both mixtures. The higher propensity of the subsonic ignition regime at lower values of $\mathit{\epsilon }$ is a result of the large IDT gradient observed in Figure 7c at T$<$T${ }_{\mathit{ave}}$.

As mentioned above, due to the absence of strong NTC effect, no coolspot is observed in the present case. Moreover, off-centered ignition only appeared for the PRF mixture in a limited zone at high values of ($\mathit{\xi }$,$\mathit{\epsilon }$), indicated by the purple dash-dotted line in Figure 7a. As apparent in Figure 7a, for the present case, the off-centered ignition region is majorly located inside the subsonic ignition regime. Moreover, the occurrence of off-centered ignition created an irregular inflection in the upper boundary of the detonation peninsula, ${\xi }_u$. Such a behavior can be justified with the increased ignition kernel size and the reduced effective distance for detonation development in the off-centered ignition mode.

Regarding the PRF-E detonation, Figure 7b shows that the detonation zone is increased compared to case III and IV. It is notable that although the detonation peninsula is located in a similar $\mathit{\xi }$,$\mathit{\epsilon }$ range for both fuels, the $\mathrm{\nabla T}$ range for PRF is double the corresponding range for PRF-E. Moreover, $r_0$ is extended to larger values for the PRF-E mixture compared to the PRF mixture.

It is observable from Figure 7(a,b) that the expected pressure distribution versus the hotspot size follows a similar pattern in the previous cases for PRF-E while it is affected by NTC chemistry for the PRF mixture. The occurrence of off-centered ignition results in a non-monotonous pressure distribution versus the hotspot radius inside the detonation peninsula for the PRF mixture. Moreover, the hotspots with similar $\mathrm{\nabla T}$ and $r_0$ for the PRF and PRF-E mixtures result in identical $P_h$ values.

Further analysis of maximum pressure distribution

In this section, the detonation zone is investigated by the second method described in Section 4.1.1, i.e. using the maximum hotspot pressure ($P_{\mathit{max}}$). In the previous section, the traveled distance by a detonation wave was an effective parameter in the recognition of detonation peninsula as well as the reported hotspot pressure $P_h$. The second approach, however, considers the maximum pressure throughout the combustion process inside the hotspot. Furthermore, the hotspots with a maximum pressure higher than the Chapman Jouguet pressure, $P_{\mathit{cj}}$, are selected as detonation zone on the ($\mathit{\xi }$,$\mathit{ϵ}$) map. Here, the pressure values are normalized by the corresponding von Neumann pressure (P_znd) for case I-V in Table 2.

In the present method, the choice of $P_{\mathit{max}}$ enables the comparison with $P_{\mathit{znd}}$. It is notable that, the hotspot-induced detonation does not possess the same structure as the ZND detonation. In the ZND detonation structure, reactions only happen behind the von Neumann spike. However, in the hotspots, detonation is initiated due to shock wave amplification by coherent energy release and therefore, reactions are partially progressed in the mixture ahead of the detonation’s leading shock front (Pan et al. 2019). In addition, the hotspot-induced deflagration to detonation followed by a large pressure intensities ($P_{\mathit{max}}$ > $P_{\mathit{znd}}$) has been reported as the mechanism of initiating super-knock in SI-engines (Wang et al. 2015). Overall, a comparison between the $P_{\mathit{znd}}$ and $P_{\mathit{max}}$ can lead to a general understanding of the strength of the hotspot-induced detonation.

Contours of $P_{\mathit{max}}$ distribution for case I and II are shown in Figure 8. The regions attributed to supersonic and subsonic ignition front propagation modes are outlined with letters P and B, respectively. A comparison between Figures 8, 3 and 4(a,b) shows that the predicted detonation peninsula is located in the similar ranges of $\mathit{\xi }$ and $\mathit{\epsilon }$ for both the PRF and PRF-E mixtures at T = 1000 K and 830 K, respectively. However, despite the positive correlation between $\mathit{\epsilon }$ and $P_h$ observed in the previous section, Figure 8 shows that $P_{\mathit{max}}$ does not necessarily increase alongside the $\mathit{\epsilon }$ axis. Once the traveled distance by the ignition front (i.e. the hotspot size) is long enough for the detonation wave development, the propagation speed converges to values close to $U_{\mathit{cj}}$ and the maximum pressure remains constant.

Figure 8. Predicted detonation peninsula of case I and II using the second method, i.e. The contours depict the maximum hotspot pressure normalized by the von Neumann pressure. Letters P, and B indicate supersonic ignition and subusonic ignition regimes.

Figure 8 shows that for all cases, $P_{\mathit{max}}$ values inside the detonation zone reach higher values than $P_{\mathit{znd}}$. This observation indicates that under certain conditions, the reaction front may form over-driven detonation inside the hotspot. Considering a planar detonation, the reason for the onset of the over-driven detonation can be an unstable (pulsating) detonation (Leung, Radulescu, and Sharpe 2010) or the onset of deflagration (or subsonic ignition) to detonation transition (DDT). As discussed in Appendix B, MSHR and the appearance of two distinct heat release pulses in the detonation profile can result in high-frequency instability of the pressure wave leading to an over-driven detonation. Nevertheless, it is observed from Figure 8 that the regions with higher $P_{\mathit{max}}$ are located closer to the border between the detonation and the subsonic-ignition zones, i.e. near the ${\xi }_u$. Therefore, a major reason for the occurrence of over-driven detonation could be the appearance of DDT as a result of the combination of the two regimes.

Figure 9. Predicted detonation peninsula of case III-V using the second method, i.E. The contours depict the maximum hotspot pressure normalized by the von Neumann pressure. Letters P, and B indicate supersonic ignition and subusonic ignition regimes.

Figure 9 shows the $P_{\mathit{max}}$ distribution inside the detonation peninsula for the high-pressure cases (case III-V). Similar to Figure 8, in Figure 9 the predicted detonation peninsula conforms with the borderlines outlined by the first method for case III-V. Moreover, Figure 9 depicts that the maximum pressure can reach higher values than the corresponding $P_{\mathit{znd}}$ values. For both fuel mixtures, it is observable that $P_{\mathit{max}}$ obtains higher values close to the upper branch of the detonation peninsula, ${\xi }_u$, indicating the DDT effect on elevating the pressure levels. However, for the PRF mixture under low and intermediate temperatures, i.e. in case IV and V, the pressure does not follow a monotonous trend. As mentioned in the previous section, NTC has a notable effect on the ignitability of the mixture, which can result in enhanced or weakened detonation inside the hotspots and consequently, higher or lower $P_{\mathit{max}}$ values.

Estimation of knock propensity

In the previous sections, the occurrence of different ignition modes inside a hotspot was analyzed. The main focus was on the detonation mode under SI-engine relevant condition. In the present section, we utilize the results of the generated regime diagrams as well as the reported results in the literature to address the relevance of the observed events to knock and super-knock in SI engines.

According to the literature, the maximum permitted knock intensity ($\mathrm{\Delta P}$) is estimated to be five bar for SI engines (Wang et al. 2015) and heavy knock has been attributed to $\mathrm{\Delta P}$ $\approx$ 24 bar (Kalghatgi and Bradley 2012). Super-knock event is designated to take place when $\mathrm{\Delta P}$ >50 bar (Pan et al. 2021; Wang et al. 2015). An additional criterion for super-knock, presented by Luong et al. (2021), classifies this event as a noticeable amount of heat release ($\stackrel{>}{\approx \phantom{{ }_2}}$ 10%) that occurs at pressure levels higher than the ZND pressure. Kalghatgi et al. (2009) used a realistic estimation of $r_0$ = 5 mm and $\mathrm{\nabla T}$ = 2 K/mm for the hotspots based on the clearance height at engines top dead center and integral length scales of turbulence and temperature gradients in promoting knock in SI engines, which is the nominal point marked with letter N in Figures 3–7.

In this section, $\mathrm{\Delta P}$(= $P_h-P_{\mathit{cv}}$) values are estimated using P_h values found in Section 4.1 on the regime diagrams, Figures 3–7 (a,b), and $P_{\mathit{cv}}$ from Table 1. Considering the pressure intensity thresholds, data from Tables 1 and 2 suggests that under low initial pressures i.e. at P = 20 bar, the occurrence of detonation may result in heavy or super-knock events since ${\mathrm{P}}_{\mathit{cj}}-P_{\mathit{cv}}\approx$ 50 bar and ${\mathrm{P}}_{\mathit{znd}}-P_{\mathit{cv}}>$ 50 bar. However, at high initial pressures, P = 50 bar, detonation can only be classified as super-knock since the difference between the detonation pressure ($P_{\mathit{znd}}$ and $P_{\mathit{cj}}$) and $P_{\mathit{cv}}$ is substantially increasing, i.e. $\mathrm{\Delta P}$ >150 bar.

Table 3 presents the calculated $\mathrm{\Delta P}$ values for detonation regime ($\mathrm{\Delta }{P}_D$), and also for the nominal point ($\mathrm{\Delta }{P}_N$) in case I-V. Next, these values are compared with the $\mathrm{\Delta }\mathit{P}$ thresholds mentioned above (Wang et al. 2015). Moreover, the predicted maximum pressure values, $P_{\mathit{max}}$, in Section 4.2 (Figures 8 and 9) are analyzed with respect to the super-knock threshold, i.e. $P_{\mathit{max}}/{\mathrm{P}}_{\mathit{znd}}$>1 according to (Luong et al. 2021).

Table 3. Pressure intensity levels at the nominal point location and inside the detonation peninsula for case I-V. Additionally, the table shows the degree of the NTC behavior observed in 1D simulations for case I-V.

	Fuel	P[bar]	Tave[K]	$\mathrm{\Delta }{P}_N$ [bar]	$\mathrm{\Delta }{P}_D$ [bar]	NTC effect
Case I	PRF	20	1000	<0	20-80	No
	PRF-E	20	1000	<0	20-100	No
Case II	PRF	20	830	30	10-70	Mild
	PRF-E	20	830	30	20-90	Mild
Case III	PRF	50	1000	90	90-190	No
	PRF-E	50	1000	<0	90-260	No
Case IV	PRF	50	900	<0	40-160	High
	PRF-E	50	900	120	60-350	No
Case V	PRF	50	830	80	50-236	Mild
	PRF-E	50	830	120	100-350	No

In case I, the nominal point lies in the subsonic ignition propagation mode for both fuels, Figure 3(a,b), therefore $P_h<P_{\mathit{cv}}$ and $\mathrm{\Delta }{P}_N<0$. Such a result indicates that hotspots with nominal conditions may not generate extra pressure fluctuations in SI engines at relatively high temperatures and low pressures. On the other hand, the values of $\mathrm{\Delta }{P}_D$, Figures 3(a,b), indicate the occurrence of heavy knock or super-knock inside SI engines where PRF-E produces stronger pressure intensities.

In case II, the nominal point is located at the edge of the detonation peninsula in Figure 4(a,b). Nevertheless, $\mathrm{\Delta }{P}_N$ values indicate that under low initial pressure and temperature conditions, hotspot ignition may result in heavy knocking for both mixtures in SI engines. It should be noted that under slightly higher temperature stratification than the nominal point, e.g. $\mathrm{\nabla }\mathit{T}$ = 3 K/mm, utilizing the ethanol blend results in a considerable reduction of $\mathrm{\Delta P}$ and yet a better performance by entering the subsonic ignition mode. Moreover, in case II, $\mathrm{\Delta }{P}_D$ values represent heavy knock or super-knock event for both mixtures.

In case III, the nominal point is located inside the detonation peninsula for PRF and in the subsonic ignition zone for the PRF-E mixture, Figure 5(a,b). Moreover, $\mathrm{\Delta }{P}_N$ values indicate that increasing the initial pressure at high temperatures may result in the occurrence of super-knock for the PRF mixture. However, PRF-E performs similar to the low initial pressure cases regarding the ignition mode; i.e. no extra pressure fluctuations in SI engines. Furthermore, $\mathrm{\Delta }{P}_D$ values inside the detonation peninsula resemble super-knock pressure intensity levels for both fuels.

In case IV, although PRF shows a broad detonation peninsula due to NTC, the nominal point is located in the supersonic ignition regime for PRF and inside the detonation peninsula for PRF-E, Figure 6(a,b). According to Table 3, the $\mathrm{\Delta }{P}_N$ values indicate low knock levels for PRF while for the PRF-E mixture the values resemble super-knock pressure intensity levels. It should be noted that an opposite trend can be observed at slightly higher temperature stratification levels, e.g. at $\mathrm{\nabla }\mathit{T}$ = 3 K/mm for both mixtures. Moreover, $\mathrm{\Delta }{P}_D$ levels indicate the super-knock event for both fuel mixtures. According to above discussions, normal knock pressure intensities may appear for the PRF mixture in the presence of supersonic ignition regime.

Finally, in case V, for both fuel mixtures the nominal point is located inside the detonation peninsula, Figure 7(a,b). The calculated $\mathrm{\Delta }{P}_N$ values indicate that the pressure intensity might correspond to super-knock. Note that for the same hotspot size, $\mathrm{\Delta }\mathit{P}$ decreases rapidly for both mixtures by increasing temperature stratification, e.g. to 3 K/mm, see Figure 7(a,b). Moreover, $\mathrm{\Delta }{P}_D$ values, inside the detonation peninsula, indicate heavy and super-knock pressure intensity levels for PRF and super-knock for PRF-E mixture.

Furthermore, considering super-knock classification by Luong et al. (2021), i.e. $P_{\mathit{max}}/{\mathrm{P}}_{\mathit{znd}}$>1, Figure 8 shows that this region is only limited to the regions at high $\mathit{\xi }$-$\mathit{\epsilon }$ for lower initial pressure cases (case I and II). However, in case II, it exhibits a vertical and horizontal shape for the PRF and PRF-E mixtures, respectively. The horizontal distribution of the super-knock region is related to larger hotspot sizes of the PRF-E mixture. The vertical distribution indicates that PRF is prone to super-knock at a wider range of $\mathrm{\nabla T}$. With respect to the higher initial pressure cases (case III-V), Figure 9 shows that in case III super-knock may be limited to high $\mathit{\xi }$-$\mathit{\epsilon }$ regions for PRF while it is extended to smaller $\mathit{\epsilon }$ values at high $\mathit{\xi }$ for the PRF-E mixture. Although in case III, considering the $\mathrm{\Delta }{P}_N$ values, the PRF-E mixture is less prone to detonation, it might produce a severe knock in the presence of detonative combustion. Moreover, Figure 9 shows that in case IV, super-knock is attributed to limited zones mostly located near the subsonic regime for both fuel mixtures. However, it is notable that the regions with $\mathrm{P}/{\mathrm{P}}_{\mathit{znd}}\stackrel{>}{\approx \phantom{{ }_2}}1$ include a relatively larger range of $\mathit{\xi }$-$\mathit{\epsilon }$ for PRF compared to the PRF-E mixture, indicating higher knock propensity of the PRF mixture. Finally for case V, Figure 9, the super-knock classified regions for the PRF-E mixture comprise a larger part of the detonation peninsula on the ($\mathit{\xi }$,$\mathit{\epsilon }$) map. However, it is notable that the associated range of $\mathrm{\nabla T}$ for PRF-E is half the range for the PRF mixture, Figure 7(a,b).

According to the above discussions:

1. Increasing the initial pressure results in strong pressure intensities (knocking) in the presence of detonation. Therefore, under-elevated pressures detonation may result in super-knock, which supports prior observations from the literature (Kalghatgi and Bradley 2012; Wang et al. 2015). Meanwhile, under-elevated initial pressures, the mechanism corresponding to the normal knock event can be relatively weak shock waves generated by the supersonic ignition regime.

2. In general, NTC chemistry increases the knock propensity by expanding $\mathrm{\nabla T}$ levels where detonation appears. However, NTC effects may suppress the pressure amplification process by reducing the detonation traveling distance while simultaneously promoting the supersonic/spontaneous ignition modes.

3. Under relatively high $\mathrm{\nabla T}$ values and low-intermediate average temperatures, PRF-E performs better than PRF under both pressure levels regarding knock resistance. However, at lower $\mathrm{\nabla T}$ levels, PRF-E may produce stronger pressure intensities.

4. Despite more limited knocking conditions and lower energy density of the PRF-E mixture, it may produce stronger pressure intensities at the onset of detonation. The present results suggest that absence of NTC chemistry in PRF-E as well as larger hotspot sizes can provide required run-up distance for developing detonation. In addition, higher deflagration tendency of the PRF-E mixture can increase the chance DDT occurrence and consequently, generation of over-driven detonation.

Analysis of the transient waves propagation

In this section, transient temperature and normalized pressure profiles ($\mathit{P}/{\mathrm{P}}_{\mathit{cj}}$) are demonstrated for selected points from the regime diagram of case IV for the PRF mixture. The intent of this section is to show the different ignition regimes and the propagation patterns that are outlined in the previous sections. We note that since the propagation patterns are similar between the different cases, here, only case IV is considered as a showcase.

In Section 4.1.5, it was discussed that NTC chemistry results in the appearance of coolspots, off-centered ignition and secondary ignition for the PRF mixture at $T_{\mathit{ave}}$ = 900 K and P = 50 bar. Figure 10 shows the location of the selected points on the ignition regime diagram of the PRF mixture in case IV. Points 1–3 are located in the supersonic ignition, subsonic ignition, and developing detonation regimes, respectively. It is notable that points 1,3 and 4 are initialized as coolspots while point 2 is a hotspot. Point 4 is located on the border between the coolspots and off-centered ignition regions and it is inside the detonation peninsula. Finally, point 5 demonstrates an example of the scenario in which NTC chemistry leads to the appearance of a secondary ignition kernel, which will be further demonstrated in the following. It is noted that the instance when the maximum temperature reaches 1800 K is labeled as time zero in the following analysis. Moreover, according to Tables 1 and 2, $P_{\mathit{cv}}/{\mathrm{P}}_{\mathit{cj}}\approx$ 0.54 and ${\mathrm{P}}_{\mathit{znd}}/{\mathrm{P}}_{\mathit{cj}}\approx$ 1.8 in case IV.

Figure 10. Locations of selected points on PRF mixture regime diagram for transient analysis at 50 bar and 900 K, case IV.

Transient reaction front analysis: point 1

Figure 11 depicts the transient temperature and pressure profiles for point 1 located in the supersonic ignition regime, initialized as a coolspot. The high propagation speeds of temperature and pressure waves are evident in Figure 11 while the estimated propagation speed is 12,620 m/s. Figure 11 shows that temperature is relatively high inside the coolspot at the instance of ignition. Moreover, the temperature difference between the fully burnt and unburnt gas inside the coolspot at each time instance is relatively low ($\mathrm{\Delta T}\approx 1000\mathit{K}$). Such characteristics have been used to distinguish between the ignition fronts and deflagration in previous studies (Karimkashi et al. 2020; Sankaran et al. 2005). It is observed in Figure 11 that the supersonic ignition results in the generation of a pressure wave of the order of $P_{\mathit{cv}}$ while approximately half of $P_{\mathit{cj}}$, see Table 2.

Figure 11. Transient profiles of temperature and normalized pressure at point 1 for PRF mixture at case IV (T = 900 K and P = 50 bar), demonstrating supersonic ignition propagation.

Transient reaction front analysis: point 2

Figure 12 shows the transient profiles of temperature and normalized pressure for a subsonic ignition front propagation (point 2). The time period from the ignition onset until when the front reaches the hotspot border shows an order of magnitude increase compared to the previous case, i.e. point 1. The estimated value for the propagation speed is 378.6 m/s, which is still much higher than the corresponding laminar flame speed of the PRF mixture ($\approx 1.3\mathit{m}/\mathit{s}$).

Figure 12. Transient profiles of temperature and normalized pressure at point 2 for PRF mixture at case IV (T = 900 K and P = 50 bar), demonstrating subsonic ignition propagation.

It is notable that the temperature difference between the fully burnt and unburnt gas inside the hotspot at each time instance in the present case is similar to that of the corresponding laminar premixed flame ($\mathrm{\Delta T}\approx 1800$ K). Moreover, the pressure profiles show that the maximum pressure is lower than both $P_{\mathit{cj}}$ and $P_{\mathit{cv}}$ throughout the whole process. Shortly after the primary ignition ($\approx$ 12 $\mathit{\mu }$s), the reactions dissociate into three separate fronts. A secondary reaction front is formed behind the main front (with a major temperature increase), generated by continuous reactions of the unburnt intermediate species that increases the peak temperature to 3000 K. Figure 12 shows that the secondary reaction front slightly lags behind the main front. Also, the pressure profiles in Figure 12 indicate the existence of a pressure wave ahead of the main reaction front. This weak pressure wave is generated by the LTI front that leads to the third inflection on the temperature profiles marked by diamonds in Figure 12. A comparison between the pressure and temperature profiles shows that the LTI front has a higher propagation speed than the main reaction front.

Transient reaction front analysis: point 3

Figure 13 shows the transient profiles of temperature and normalized pressure for a detonative case (point 3), initialized as a coolspot. As it can be observed from Figure 10, the simulation point on the map is located slightly outside of the detonation peninsula and close to the lower detonation zone boundary. However, the temperature and pressure profiles in Figure 13 indicate the formation of a detonation front. Particularly, Figure 13 depicts that the front pressure finally reaches the CJ pressure level inside the coolspot. However, the pressure magnitude is still lower than the von-Neumann pressure. As it can be seen from the profiles, the propagation speed is higher during the initial stages. The low reactivity distribution inside the coolspot, due to NTC chemistry, leads to supersonic/spontaneous ignition front propagation during the early stages. Therefore, the ignition front does not turn into a fully developed detonation inside the coolspot, which results in a deficit between the front and the $P_{\mathit{znd}}$ pressures. Figure 10, shows that increasing the coolspot size leads to a higher $P_h$ value and the corresponding point would then locate inside the detonation peninsula.

Figure 13. Transient profiles of temperature and normalized pressure at point 3 for PRF mixture at case IV (T = 900 K and P = 50 bar), demonstrating detonation development.

Transient reaction front analysis: point 4

Figure 14 depicts the temperature and normalized pressure evolution for an off-centered detonation case at point 4. As described in Figure 10, point 4 is initialized as a coolspot and the onset of the ignition occurs with an offset from the center of the coolspot. From Figure 14 it is observed that the left and right propagating ignition fronts evolve into detonation waves. A comparison between points 3 and 4 indicates that they have the same coolspot size while the coupling factor $\mathit{\xi }$ is higher for point 4. Subsequently, compared to point 3, it is observed that the transition from the initial supersonic/spontaneous propagation to detonation happens in a more prompt manner. As a result, the left propagating detonation wave can reach the fully developed stage with pressure levels above the CJ pressure. Moreover, from Table 2, ${\mathrm{P}}_{\mathit{znd}}/{\mathrm{P}}_{\mathit{cj}}\approx 1.8$, which indicates that pressure reaches levels higher than the ZND detonation pressure in the present case. In addition, the profiles of the left propagating detonation depict variations in values lower and higher than ${\mathrm{P}}_{\mathit{znd}}/{\mathrm{P}}_{\mathit{cj}}$, indicating an unsteady (pulsating) detonation front. Such pressure variations could be originated from the MSHR characteristics of the fuel mixture detonation as depicted and discussed in Figure B1.

Figure 14. Transient profiles of temperature and normalized pressure at point 4 for PRF mixture at case IV (T = 900 K and P = 50 bar), demonstrating occurrence of off-centered ignition.

Transient reaction front analysis: point 5

Figure 15 shows the transient profiles of temperature and normalized pressure for a detonative case at point 5 in Figure 10. The temperature profiles depict that the primary ignition at the center forms a subsonic ignition front as well as an LTI wave ahead of it propagating toward the right end of the domain. Similar to point 2, the LTI front increases temperature and pressure and at t = 15.4$\mathit{\mu s}$, it forms a secondary ignition kernel. Later, the secondary ignition kernel rapidly develops into right and left propagating detonation waves, as depicted in Figure 15. As mentioned earlier in Section 4.1.5, NTC chemistry can change the dynamics of reaction waves and result in a complex system of shock-reaction inside the hotspot. Moreover, the pressure profiles indicate that the generated detonation wave from the secondary ignition kernel increases the pressure to levels higher than $P_{\mathit{cj}}$.

Figure 15. Transient profiles of temperature and normalized pressure at point 5 for PRF mixture at case IV (T = 900 K and P = 50 bar), demonstrating appearance of secondary ignition.

As a final remark, the current work presents a detailed study on abnormal combustion modes inside to SI engines conditions through exploring various possible ignition regimes using two different gasoline surrogates. In this work, ignition regime diagrams were presented within a similar framework of the original detonation peninsula (Bradley et al. 2002) for the mixtures with NTC effect at different average temperatures, pressures, temperatures stratification levels and radii. In addition, estimation of the hotspot pressure levels provides the required information to relate the present work results to knock/super-knock in SI engines using available definitions from the literature. We note that the utilized numerical and theoretical framework (Bradley et al. 2002) are particularly designed for 1D spherical hotspot simulations. However, by extension of the numerical simulations to 2d and 3d, the impact of physical effects such as turbulence, shock reflection, and cellular structures of detonation could be modeled accordingly using the present reactive flow solution approach. As noted in Section 2.2.3, the reflection of the auto-ignition induced shock waves from the right boundary is not considered in present work. Such an event can interact with the system of shock/reaction fronts inside the hotspot (Dai et al. 2021; Terashima and Koshi 2015). Nevertheless, such simulations are left for future work and they are not further considered herein.

Conclusions

The present work continues the authors’ previous study on the effect of thermal stratification on the prevalence of deflagration and ignition modes under SI engine relevant conditions. Here, we study the different ignition regimes and their impact on hotspot pressure levels in the end-gas region. The work utilizes a 1D numerical and theoretical framework to generate ignition regime diagrams (detonation peninsulas) based on reactivity ($\mathit{\epsilon }$) and coupling factor ($\mathit{\xi }$) parameters. The regime diagrams are generated for different initial pressures (20 and 50 bar) and average temperatures (830–1000 K) for sensitive (PRF-E) and nonsensitive (PRF87) gasoline surrogates. Furthermore, available data on hotspot pressure levels are related to different knocking scenarios in SI engines using available definitions from the literature. The following conclusions are made:

1. The generated detonation peninsulas confirm the sensitivity of the ignition regimes to the choice of fuel mixture, initial pressure, and average temperature. Moreover, the present work denotes that multistage heat release and NTC behavior of gasoline surrogates can change the distribution of the ignition regimes inside the ($\mathit{\xi }$,$\mathit{\epsilon }$) maps.

2. NTC may have both prohibiting or inhibiting impact on detonation development. NTC chemistry is noted to extend the supersonic/spontaneous ignition as well as detonation regimes to higher temperature stratification levels. On the other hand, NTC may suppress the detonation development by reducing the traveling distance of the reaction front, by initiating secondary ignition kernels and off-centered ignition.

3. A limited ignition and detonation propensity was observed for the PRF-E mixture compared to PRF due to lower NTC behavior. In other words, ethanol addition to gasoline is noted to reduce the knock-propensity of the mixture. The lower energy density and longer excitation times mitigate the generated pressure perturbations and thus reduce the formation of intense shock waves for the PRF-E mixture. However, the probability of heavy knock and super-knock increases at the onset of detonation inside the hotspots of the PRF-E mixture. Longer excitation times and a higher tendency for a blended deflagration-detonation mode for PRF-E may provide required run-up distances for detonation development and thus, increasing the probability of DDT.

4. Ignition delay time (${\tau }_i$) has the highest impact on hotspot pressure levels of the studied gasoline surrogates, compared to energy density ($E_c$) and excitation times (${\tau }_e$). First, ${\tau }_i$ has a direct impact on the dynamics of ignition front propagation, ignition kernel size, and the traveled distance by the ignition front. Second, there is a strong correlation between ${\tau }_i$ gradient ($\mathit{\alpha }$) and $\mathit{\xi }$.

5. Under low initial pressures, engine heavy knocking may be connected to the onset of detonation. On the other hand, for high initial pressures, heavy knock can be related to supersonic ignition propagation. However, regardless of the fuel type (octane sensitivity), at higher pressures, onset of detonation may be attributed to the super-knock event.

Acknowledgement

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 [grant numbers 318024 and 332784].

Appendix A.

Solver validation

In this section, we present the numerical solver validation. Initially the capability of the solver in preserving the discontinuities in non-reacting compressible flows in a shock tube is examined. Furthermore, as a reacting detonative case, hotspot detonation of methane is simulated and results are compared with available simulation results from Pan et al. (2019). Finally, grid size independency is assessed for the main hotspot-induced ignition simulations.

A1. Non-reacting shocktube

For the non-reacting test, sod shocktube (Sod 1978) problem is considered to evaluate the solver’s capability in predicting the displacement of discontinuities. The sod shocktube consists of a one dimensional Riemann problem. Therefore, the simulation results can be compared with the analytical solution of the problem. This makes sod’s problem an ideal case for validation of the numerical methods. The simulation case is a semi-infinite 1D shock tube, where the left and right boundaries are considered to be far away enough to avoid reflections at the ends. The gas inside the shock tube is a pure inviscid gas with a molecular weight of 28.9, and it is separated initially by a diaphragm in the middle. Flow proprieties, i.e. velocity (U), pressure (P) and density $\mathit{\rho }$, at the left and right sides of the diaphragm are mentioned in Table A1.

Table A1. Test case configuration.

${\rho }_{\mathit{l}}$ [Kg/m^3]	${\rho }_{\mathit{r}}$ [Kg/m^3]	Pl [pa]	Pr [pa]	U [m/s]	MW	L [m]	Gridpoints
1	0.125	1	0.1	0	28.9	1	200-800

The initial discontinuity later transforms into a right propagating shock wave, a left propagating rarefaction wave and a contact discontinuity moving to the right side with a lower speed than that of the shock wave. Figure a1 depicts pressure, velocity, density and acoustic speed distributions for the reactor compared against the analytical solution at 0.1644 s for different grid resolutions. The convective courant number is equal to 0.02 for these simulations (the acoustic courant number $\approx$ 0.04).

Figure A1. Sod’s shocktube pressure, velocity, density and acoustic speed distribution versus analytical solution at 0.1644 s.

As it can be seen from Figure A1, the simulation results compared to the analytical solution indicate that the discontinuities are well captured with the solver with various grid resolutions.

A2. CH4 Hotspot detonation simulation

The second case investigates the performance of the solver in predicting the detonation inside a hotspot, modeled as a localized 1D temperature gradient inside a constant volume reactor. The simulation setup is similar to the unsteady simulations mentioned in the paper. Results from numerical simulation on methane hotspot detonation by Pan et al. (2019) is considered for the validation purpose. The simulation setup presented in Table A2, i.e. domain length (L), grid size ($\mathrm{\Delta x}$), hotspot radius ($r_{\mathit{hs}}$), average temperature ($T_{\mathit{ave}}$), temperature gradient ($\mathrm{\nabla T}$), initial pressure ($P_0$) and $\mathit{\xi }$, is identical to those in Ref (Pan et al. 2019).

Table A2. Test case configuration.

Tave [K]	$\mathrm{\nabla }\mathit{T}$ [K/m]	P0 [atm]	rhs [m]	$\mathit{\xi }$	L [m]	$\mathrm{\Delta x}$ [m]	UCJ [m/s]	HRL [m]
1381	1000	40	0.008	1	0.04	8e^−7	1834	3e^−5

The GRI 3.0 (Smith 1999) detailed mechanism has been considered as the kinetic model for CFD simulations, as well as for the calculation of $U_{\mathit{cj}}$ and HRL, Table A2, using the SDToolbox (Kao and Shepherd 2008). It is worth noting that due to the sensitivity of $\mathit{\xi }$ values to the choice of the mechanism (Pan et al. 2021; Su, Dai, and Chen 2021), this value is reported from Ref (Pan et al. 2019).

Figure A2 shows the temperature profiles at 136–139 $\mathit{\mu }$s inside the hotspot as well as the same profiles from Ref. (Pan et al. 2019). The transient temperature profiles are in agreement with the reference data. From Table A2, it is observable that detonation structure is well resolved within the chosen mesh with $\approx$ 37 cells per HRL. Moreover, the propagation speed estimated by tracking the pressure and reaction fronts is plotted compared with the calculated CJ speed. It is observable that shortly after the initial ignition (at $\mathit{r}/r_0\approx 0.3$), the pressure and reaction waves merge and form a detonation front with a propagation speed close to the CJ speed.

Figure A2. Temperature profiles (dashed lines are from (Pan et al. 2019), solid lines are from the present work) and reaction, pressure fronts speeds for methane hotspot detonation, $P_0=40$ atm, $T_{\mathit{ave}}=1381$ K, $\mathrm{\nabla T}=1$ K/mm, $r_0=8$ mm. CJ speed is calculated at $T_{\mathit{ave}}$ using SDToolbox (Kao and Shepherd 2008).

A3. Grid size sensitivity analysis

In the present work, a uniformly refined mesh is utilized inside the hotspot for resolving the evolution of ignition fronts inside the hotspots. The choice of suitable grid sizes is based on sufficient resolution of the zone between the leading shock front and the peak of heat release, $L_2$L in Figure B1, in ZND detonation structure. Based on the depicted ZND profiles in Figure B1, for cases studied in this work the minimum characteristic length is attributed to case III, i.e. the PRF mixture at P = 50 bar and $T_{\mathit{ave}}$ = 1000 K where $L_2\approx 1e^{-5}\mathit{m}$. Table A3 demonstrates the specific hotspot configuration considered for the grid size sensitivity analysis, i.e. domain length (L), grid size ($\mathrm{\Delta x}$), hotspot radius ($r_0$), average temperature ($T_{\mathit{ave}}$), temperature gradient ($\mathrm{\nabla T}$), initial pressure ($r_0$), $\mathit{\xi }$ and $\mathit{ϵ}$.

Table A3. Grid sensitivity case configuration.

Tave [K]	$\mathrm{\nabla }\mathit{T}$ [K/m]	P0 [bar]	R0 [m]	$\mathit{\xi }$	$\mathit{\epsilon }$	L [m]	$\mathrm{\Delta x}$ [$\mathit{\mu }$m]	HRL [m]
1000	502	50	0.004	4	10	0.04	0.5, 0.75, 1	1e^−5

Figure A3 depicts the normalized pressure (P/Pznd), temperature and OH mass fraction transient profiles. As it is observable from Figure A3, at the mesh size of 1 $\mathit{\mu m}$ (the black dashed-line), the waves travel slightly further distances compared to 0.5 $\mathit{\mu m}$, and 0.75 $\mathit{\mu m}$ mesh results. Moreover, the peak pressure at the last time step is under-predicted. However, the profiles for 0.5 and 0.75 $\mathit{\mu m}$ grid resolutions conform well with each other which indicates the convergence of the results with respect to the chosen grid size.

Figure A3. Mesh sensitivity study for case.Black dashed-line result from grid size = $1\text{mum}$, Blue solid line results from grid size = 0.75 $\text{mum}$, red dashed-line results from grid size = 0.5$\text{mum}$.

Appendix B.

ZND planar detonation profiles

Figure B1. The ZND detonation structures of case I-V for PRF and PRF-E mixtures. The blue lines depict the heat release profiles, and the red-dashed lines depict the temperature profiles. L1 and L2 indicate the first and second induction lengths, respectively.

For all the studied cases, Figure B1 shows a two-stage heat release profile as well as two inflection points in the temperature profile. The first inflection point is located after the initial decrease in temperature due to the endothermic C-H bond-breaking reactions. Furthermore, the major increase in temperature occurs after the second inflection point, i.e., the location of the second heat-release peak, Figure B1.

Appearance of multi-stage heat release (MSHR) has been reported in previous experiments and in ZND simulations of some fuel blends with NO₂, such as H₂/NO₂/N₂O₄ (Davidenko, M´evel, and Dupr´e 2011; Joubert, Desbordes, and Presles 2008) and CH₃NO₂/O₂ (Sturtzer et al. 2005) and CH₃NO₂/O₂ sturtzer2005origin. Moreover, such a phenomenon has been observed for the DME/O2 (Ng et al. 2009) mixture at both high (M´evel and Gallier 2018; Ng et al. 2009) and low post-shock temperatures (Han et al. 2021). Furthermore, the MSHR has been observed for large hydrocarbons with LTC (Griffiths and Scott 1987; Han et al. 2021) as well as in ZND detonation profiles of the $\mathit{n}$-heptane/CO₂/O₂ mixture at a low post-shock temperature (Liang, M´evel, and Law 2018). The two-stage heat release has been attributed to”double cellular” detonation structure, where secondary cells of smaller sizes appear inside the primary detonation cells (Ng et al. 2009; Sturtzer et al. 2005). These studies reported an order of magnitude longer induction length for the second peak of heat-release ($L_2$) than the first one ($L_1$). The presence of two distinct chemical length scales can result in high-frequency instability of the pressure wave leading to an over-driven (pressure/velocity higher than the corresponding CJ value) detonation. It should be noted that the double cellular structure has been reported to originate from the strong intrinsic instability of the unstable detonation of lean mixtures of DME/O2 as well (M´evel and Gallier 2018; Ng et al. 2009). Therefore, we acknowledge that 2D-3D simulations and detailed chemistry analyses are required to confirm the existence of the double cellular structure for gasoline surrogates under the studied conditions.

The ratio between the first and second induction length ($L_1/L_2$L) in the present cases is reported in Table 2. It is notable that for the studied cases, $L_2$ (also introduced as HRL in Section 2.2.3) is at least an order of magnitude longer than $L_1$, c.f., Figure B1 and Table 2. According to the above discussions, such a characteristic feature indicates the potential instability of detonation in both fuel mixtures. Moreover, Figure B1 depicts that at lower initial pressures, case I (II), the relative first/second stage heat release is higher than the corresponding ratio at higher initial pressures, case III (V). Such a trend denotes that the dynamics of the detonation wave propagation is majorly controlled by the second stage heat release at higher pressure levels (Ng et al. 2009).

Appendix C.

Hotspot pressure calculation

The approach for determining the detonation pressure in the first method utilized in this work is presented herein. This approach is based on clustering the maximum pressure values into a relatively high level, attributed to the detonation wave, and a low level, related to the supersonic/spontaneous and subsonic ignition modes. Using the Scikit-learn Python module (Pedregosa et al. 2011), the Mean Shift algorithm (Comaniciu and Meer 2002) is chosen for clustering the maximum pressures data. Figure C1 depicts the transient profiles of temperature and pressure as well as the maximum values of pressure profiles divided into two clusters.

Figure C1. Transient temperature and pressure profiles as well as distribution of maximum pressure for a hotspot detonation scenario. Clusters are separated by colors. Yellow circles indicate the mean value calculated for each cluster.

In Figure C1, the upper and lower clusters are colored by blue and green dots. Furthermore, the calculated mean values are indicated by yellow circles, for upper ($P_1$) and lower ($P_2$) clusters. Moreover, the attributed normalized length of each cluster is shown by $L_1$ and $L_2$ for the upper and the lower cluster, respectively.

The procedure depicted in Figure C1 is taken due to the fact that combustion progresses at different rates for the studied ignition modes and therefore, time (ensemble) average is biased by the slower process, with more data points. The issue is tackled by taking the spatial distribution of the maximum pressure into account. From Figure C1, it is observable that the lengths attributed to each cluster (the normalized radius swept by high/low pressure reaction fronts) can be then used to calculate the weighted average of the maximum pressure where $P_h=P_1\times L_1+P_2\times L_2$.