Analysis of the combustion speed in a spark ignition engine fuelled with hydrogen and gasoline blends at different air fuel ratios

Full Length Article Analysis of the combustion speed in a spark ignition engine fuelled with hydrogen and gasoline blends at different air fuel ratios Pedro Gabana a,* , Blanca Gim´ enez a , Jos´ e Martín Herreros b , Athanasios Tsolakis b a Department of Energy and Fluid Mechanics Engineering, University of Valladolid, Paseo del Cauce 59, E-47011 Valladolid, Spain b Department of Mechanical Engineering, University of Birmingham. Edgbaston, Birmingham, B15 2TT, UK ARTICLE INFO Keywords: Hydrogen-Gasoline combustion Quasi-dimensional combustion diagnosis Turbulent combustion speed Premixed combustion Combustion speed analysis Expansion speed ABSTRACT The use of hydrogen in internal combustion engines is a promising solution for the decarbonisation of the transport sector. The current transition scenario is marked by the unavailability and storage challenges of hydrogen. Dual fuel combustion of hydrogen and gasoline in current spark ignition engines is a feasible solution in the short and medium term as it can improve engine efficiency, reduce pollutant emissions and contribute significantly in tank to wheel decarbonisation without major engine modification. However, new research is needed to understand how the incorporation of hydrogen affects existing engines to effectively implement gasoline-hydrogen dual fuel option. Understanding the impact of hydrogen on the combustion process (e.g. combustion speed) will guide and optimize the operation of engines under dual fuel combustion conditions. In this work, a commercial gasoline direct injection engine has been modified to operate with gasolinehydrogen fuels. The experiments have been carried out at various air–fuel ratios ranging from stoichiometric to lean combustion conditions at constant engine speed and torque. At each one of the 14 experimental points, 200-cycle in-cylinder pressure traces were recorded and processed with a quasi-dimensional diagnostic model and a combustion speed analysis was then carried out. It has been understood that hydrogen mainly reduces the duration of the first combustion phase. Hydrogen also enables to increase air excess ratios (lean in fuel combustion) without significantly increasing combustion duration. Furthermore, a correlation is proposed to predict combustion speed as a function of the fuel and air mixture properties. This correlation can be incorporated to calculate combustion duration in predictive models of engines operating under different fuel mixtures and different geometries of the combustion chamber with pent-roof cylinder head and flat piston head. 1. Introduction The transport sector is increasingly moving towards decarbonisation. This sector is largely sustained by the use of fossil fuels that are burned in gas turbines or internal combustion engines (ICEs) to obtain mechanical energy. Within the transport sector, ICEs are the most widespread propulsion systems because they can work in a wide range of operating conditions. On the other hand, the lower the carbon content of fuels, the lower the CO2 emissions in relation to the energy produced in their combustion [1]. As hydrogen has no carbon atoms in its composition, it does not generate CO2 during combustion. Therefore, hydrogen is presented as an alternative to fossil fuels, since its use as a fuel in engines means the reduction of CO2 [2], CO, unburned hydrocarbons and soot emissions [3]. The use of synthetic fuels produced from renewable energy is a way to move away from the use of fossil fuels in ICEs. H2 can be produced from renewable electricity with zero CO2 equivalent emissions [3]. Hydrogen in the mobility sector is majorly used in two ways. The first one is the use of hydrogen in fuel cells [3,4]. The second option is the use of hydrogen in thermal engines [3,5]. Within the reciprocating thermal engines, the use of hydrogen in compression ignition engines (CIE) has the disadvantage that the autoignition temperature of hydrogen is higher than those of conventional fuels [5]. In spark ignition engines (SIE), the use of pure H2 as fuel enables the use of high air–fuel ratios and high levels of residual gas concentrations. The higher the air–fuel ratio or residual gas concentrations, the lower the in-cylinder temperature during the combustion process. The lower the combustion * Corresponding author. E-mail addresses: [email protected] (P. Gabana), [email protected] (B. Gim´ enez), [email protected] (J. Martín Herreros), [email protected] (A. Tsolakis). Contents lists available at ScienceDirect Fuel journal homepage: www.elsevier.com/locate/fuel https://doi.org/10.1016/j.fuel.2024.133563 Received 18 August 2024; Received in revised form 11 October 2024; Accepted 25 October 2024 Fuel 381 (2025) 133563 Available online 3 November 2024 0016-2361/© 2024 The Authors. Published by Elsevier Ltd. This is an open access article under the CC BY license ( http://creativecommons.org/licenses/by/4.0/ ). temperature, the lower the NOx emissions. Performance improvements have also been observed as reported by Kapus et al. [6]. However, the implementation of hydrogen presents several challenges, some of the most important ones are storage and availability [2,7]. The use of dual fuel (hydrogen and gasoline) is presented as a shortterm and medium-term alternative to reduce CO2 emissions by unit of energy produced [8,9] without major engine modifications. Previous works have studied how the addition of hydrogen to gasoline affects engine operation [9–11]. Like engines running only with hydrogen, the dual combustion of hydrogen and gasoline in SIE allows NOx emissions reductions by operating the engine at higher air–fuel ratio [12]. When hydrogen is added to gasoline, higher combustion speeds are obtained compared to those obtained when the engine operates only with gasoline [13]. It also enables the air–fuel ratio [12] and residual gas concentrations to be increased while maintaining a stable engine operation. Verhelst et al. [5] pointed out that more studies of hydrogen combustion under engine conditions were required. Until now, the lack of this type of studies is particularly noticeable for the case of hydrogengasoline spark ignition engines. The methodology used to study the combustion process characteristics in hydrogen-gasoline blends can be used in the combustion analysis of hydrogen blends with other renewable fuels and under other different engine operating conditions. Going beyond, studying combustion speed permits determine the thermo-fluidmechanic responsible on how the combustion processes take place unveiling new knowledge to understand the impact of fuels on combustion characteristics. Applying thermodynamic diagnostic models to in-cylinder pressure traces calculates the amount of energy released in a certain angular timing. The rate at which energy is released is a consequence of the product of the combustion speed and the flame front area. To calculate the combustion speed from heat release rate, it is necessary to use a geometrical model to determine flame front area [14]. The combustion speed calculated with both thermodynamic and geometrical models can be correlated with thermo-fluid-mechanical variables inside the combustion chamber (laminar combustion speed, turbulence intensity, integral length scale, flame front thickness, expansion speed and flame front perimeter). This knowledge translated into a correlation can be used to predict in-cylinder pressure from the thermo-fluid-mechanical variables in the combustion chamber. For this purpose, the combustion speed obtained with the correlation is used as input to the thermodynamic and geometric models. This work studies the effects on combustion speed of substituting part of the gasoline fuel with H2 in a commercial spark ignition engine. For this purpose, a blend of gasoline and ethanol is experimentally investigated as the reference/baseline fuel being partially replaced by H2. The in-cylinder pressure records are obtained and processed within a two-zone thermodynamic diagnostic model [14]. The results of the Nomenclature Variables AArea (m2) clTurbulence dissipation constant ( − ) cpSpecific heat (J/kg/K) ctTurbulence production constant ( − ) DiMass diffusivity of species iin the mixture (kg/m2) HpHeat power (J/kg) kturbulent flux energy (J) KMean flux energy (J) kδlFlame front thickness multiplication term ( − ) LMarkstein length (m) LLength scale (m) Leeff Effective Lewis ( − ) mMass (kg) pPressure (Pa) PTurbulence production term ( − ) PfFlame front perimeter (m) ˙ QHeat transfer (J/s) RRadius (m) SCombustion speed (m/s) TTemperature (K) uʹTurbulence intensity (m/s) VVolume (m3) Acronyms BMEP Brake Mean Effective Pressure (bar) CA Crank Angle (◦) CIE C ompression Ignition Engine CoV Coefficient of Variation EVC Exhaust valve closing EVO Exhaust valve opening FSR Flame speed ratio GGasoline H2 Hydrogen ICE Internal Combustion Engine IVC Intake valve closing Acronyms (cont.) IVO Intake valve opening MGasoline hydrogen mixture tests MFB Mass fraction burned NG Natural Gas SIE Spark Ignition Engine SoI Start of combustion Greek letters α Angle (◦) α uThermal diffusivity unburned products (m2/s) βZel’dovich number ( − ) γspecific heats ratio ( − ) δFlame front thickness (m) λAir fuel equivalence ratio ( − ) λTThermal conductivity (J/(m s K)) ρ Density (kg/m3) σ Density ratio ( − ) ω Amplitude growth rate ( − ) Subscripts bBurned dDiagnostic eExpansion fFlame front iIntegral kKolmogorov lLaminar max Maximum pPiston pu Unburned products in contact with the piston ref Reference res Residual gas concentration tTurbulent uUnburned wWoschni P. Gabana et al. Fuel 381 (2025) 133563 2 thermodynamic model are fed into a newly developed geometrical model for pent roof cylinder head combustion chamber which is defined by 6 parameters enabling to adapt the model to cylinder heads with different angle between valves, bores and intersection with the engine head gasket. The outputs of the geometrical model are used to predict combustion speed as a function of the properties of the fuel–air mixtures by using a new methodology, underpinned by Gim´ enez et al. [14] (NG and H2 blends in a Heron piston top combustion chamber). A new methodology for calculating and fitting the correlation of the combustion speed has been proposed. This work demonstrates that the newly proposed correlation can be used to predict the combustion characteristics for different air–fuel ratios and new fuels (e.g. gasoline-hydrogen blends), which could be extended to a variety of spark ignition combustion chamber geometries. 2. Methodology 2.1. Experimental setup The engine studied is a spark ignition Gasoline Direct Injection Engine. The engine specifications are presented in Table 1. The original engine has been modified to add H2 in the intake pipe of the engine. The H2 flow rate was measured with a H2 volumetric flow meter (Platon NG series VA GTF-2AHD with a range of 2–44 L/min 1.013 bar 20 ◦C). Thus, a mixture of air and hydrogen is introduced into the cylinder. The gasoline is supplied directly into the cylinder. Combustion data was measured and logged in one cylinder using an AVL piezo-electric in-cylinder pressure transducer, a charge amplifier (AVL FlexIFEM) and a Baumer (720 pulse per revolution) magnetic crank angle encoder. The pressure data have been recorded each 0.1 crank angle interval using an AVL Indicom system. A schematic diagram of the experimental setup is shown in Fig. 1. 2.2. Experimental methodology Table 2 shows an overview of the test plan for the engine operating points. Two fuels have been used, gasoline and ethanol E10 (G) and the same fuel being replaced by hydrogen at 10 % by mass and 87.3 % by volume (M). The fuel–air ratio was also modified from stoichiometric to λ=1.2 in Gand to λ=2 in M. All experiments have been carried out keeping the engine speed at 2000 rpm and the brake torque at 40 Nm. The spark timing is set in order to have a mass fraction burned of 50 % (CA50) in a crank angle (CA) of between 3 and 4◦after top dead centre as in previous works [12]. At each experimental point, 200 cycles have been recorded and processed. The fuel flow has been tuned to achieve the target effective torque (40 Nm). By increasing air–fuel equivalence ratio (λ), the intake pressure has been gradually increased to introduce the air or H2/air mixture (depending on the case) required to achieve the λvalue specified. Consequently, the maximum engine pressure is increased. The amount of energy introduced into the engine per cycle is similar in all experimental points. As λincreases, less fuel is introduced, this is because the indicated efficiency increases with λ. 2.3. Data processing The methodology followed during the diagnostic process is summarised in the following steps: 1. Data filtering using an adaptive filter based on Savitzky-Golay polynomials [15] to obtain a continuous and derivable pressure signal at each cycle. 2. Determination of pressure and angle offsets, compression ratio (10.8) and heat transfer coefficient (0.9). For this, a methodology based on genetic algorithms is used, as stated by Reyes et al. [16]. Except for the pressure offset, which is chosen for each cycle, the rest of the variables take the same value in all experimental conditions. 3. Processing of the filtered pressure data with a two-zone diagnostic thermodynamic model [14,17]. The thermodynamic diagnostic model allows calculating the burned mass rate, dmb/dt, from dp/dt obtained processing pressure data. More information about the thermodynamic model is available in Mendeley data section. 4. Processing of the resulting data of burned products volume with a geometric model that enables calculating flame front radius (Rf) and flame front area (Af). The flame front area is combined with dmb/dt to calculate the combustion speed (Sd) Eq. (1). dmb dt =SdAf ρ u(1) 5. Average combustion speed calculation for each flame front position in a population of 200 cycles. A cycle averaged combustion speed result is obtained for each test point and flame front position inside the combustion chamber. The results are independent of the crank angle at which they take place and dependent on the position in the combustion chamber. 3. Geometric model Based on the volume of mass burned, which is the result of the thermodynamic model, the idea is to determine the radius of the sphere centred on the spark plug intersecting with combustion chamber walls and the flat piston. The combustion chamber is defined by the 6 parameters. Once the radius of the sphere has been determined, the next objective of the model is to calculate the area of the sphere contained in the combustion chamber. This area is the one used in Eq. (1). For this purpose, the combustion chamber in the cylinder head is assumed to be delimited by (Fig. 2): - Two inclined and symmetric planes forming the angle α t, the intersection horizontal line of the two planes is located at a distance Zv from the spark plug. - A half-angle cone α Cwith a base at the cylinder head gasket having a diameter equal to the Bore. - Two vertical planes at a distance Ypfrom the cylinder axis. - A horizontal plane at a distance Ztover the spark plug. - The plane of the cylinder head gasket is at a distance Zcfrom the spark plug. The combustion chamber volume below the cylinder head gasket is delimited by: - Cylinder head gasket plane. - The cylinder of diameter 2 Rp. - The piston plane located at a distance Xpfrom the cylinder head gasket. Table 1 Engine geometric specifications. Bore x Stroke 84 x 90 mm Number of cylinders 3 Number of Valves 12 Displacement 1497 cm 3 Compression Ratio 11:1 IVO 20 CAD IVC −88 CAD EVO 92 CAD EVC 0 CAD Inlet valve peak lift 10 mm Exhaust valve peak lift 9 mm Total Crank/Piston offset 9.6 mm P. Gabana et al. Fuel 381 (2025) 133563 3 The volume covered by the sphere of radius Rfis calculated from horizontal planes at different distances Zfrom the spark plug. The intersections of these planes with the combustion chamber result in a geometry of stadiums defined by two parallel lines and a circle centred on the spark plug. The intersection of the plane with the sphere volume has the same geometry unless the sphere does not intersect in that plane with the cylinder head, then it is a circle. The burned products volume outside the cylinder head corresponds to the intersection of a sphere of radius Rfwith the plane of the piston and cylinder. Full geometric model description is available in Mendeley data section. In Fig. 2 two side views with the parameters characterizing the model are shown. The floor plane shows a horizontal plane located at a distance Zfrom the spark plug. Fig. 3 shows an image of the burned volume at a certain piston position and for a certain flame front radius. Fig. 3 has been obtained with Autodesk Inventor software. This software has been used to validate the geometric model. 4. Combustion properties 4.1. Until this section, the methodology for processing the experimental data to calculate the combustion speed (Sd) has been presented. These results and their analysis are presented in section 5.1. The final objective of this work, as mentioned in the introduction, is to predict combustion speeds (St) from the thermo-fluid-mechanical variables inside the combustion chamber. This section explains how to calculate the variables (laminar combustion speed, residual exhaust gases concentration, turbulence intensity and flame properties) dependent on the conditions in the combustion chamber necessary to fit the correlation of the combustion speed (St) in section 5.2.Laminar combustion speed Hydrogen-gasoline mixtures laminar combustion speed (Sl) determination has been previously addressed in the literature [18,19]. The laminar combustion speed depends on the pressure and temperature conditions and on the mixture composition, determined by the proportion of hydrogen in the fuel, the air excess ratio (λ) and mass fraction of combustion products Yres. The expression proposed by [19] has been used as Eq. (2), which is valid for the engine operating conditions. Fig. 1. Experimental facility schematic diagram. Table 2 Characteristics of the 14 test points analysed at Me =40Nm, n =2000 rpm, BMEP =3.35 bar. Name λ˙ mf(kg/h) ˙ ma(kg/h) SoI (◦CA) pIVC(bar) pmax(bar) Indicated efficiency (%) Yres(%mass) G1 1.0 0.785 11.614 –32.4 1.01 29.51 34.9 19.9 G1.1 1.1 0.760 12.369 −35.9 0.95 29.49 36.4 18.9 G1.2 1.2 0.760 13.481 −41.2 1.01 29.87 36.5 18.4 M1 1.0 0.638 10.833 −20.3 1.02 28.87 34.9 21.7 M1.1 1.1 0.657 12.283 −21.1 1.11 30.09 35.8 20.7 M1.2 1.2 0.634 12.942 −21.7 1.12 30.04 36.9 21.2 M1.3 1.3 0.623 13.784 –23.6 1.17 30.67 37.7 21.2 M1.4 1.4 0.612 14.580 −25.2 1.21 30.90 38.5 21.2 M1.5 1.5 0.602 15.389 −26.5 1.24 31.57 39.2 21.2 M1.6 1.6 0.610 16.605 −28.5 1.28 32.00 39.7 21.2 M1.7 1.7 0.593 17.172 −28.5 1.30 31.66 40.2 21.2 M1.8 1.8 0.588 18.040 –33.1 1.37 33.18 40.5 21.1 M1.9 1.9 0.585 18.945 −36.2 1.40 33.30 40.7 20.7 M2 2.0 0.586 19.984 −41.4 1.45 34.54 41.1 20.6 P. Gabana et al. Fuel 381 (2025) 133563 4 Sl=Sl0(Tu Tref ) α T(p pref )βp (1−2.1Yres)(2) Where Tu(K) is the unburned products temperature, p(bar) is the pressure, and the reference temperature and pressure are respectively 298 K and 1 bar. Sl0, α Tand βpcan be calculated according to the expression proposed in [19]. These three parameters depend on the hydrogen volume fraction in the intake gas and λ. Eq. (2) provides laminar combustion speed (Sl) dependent on the combustion chamber conditions. Experiments at low pressures and temperatures show that the combustion speeds of ethanol and iso-octane are similar [20]. Simulations for the calculation of the laminar combustion speed of gasoline-ethanol blends at high pressures and temperatures show that the presence of ethanol in gasoline has little influence on the laminar combustion speed compared to pressure, temperature, fuel air ratio and recirculated exhaust gases [21]. For this reason, the laminar combustion speed has been considered to be similar to a mixture of iso-octane and hydrogen; however, the flame properties in section 4.4 have been calculated considering the composition of the mixture. 4.2. Residual exhaust gases The residual gases mass fractions were calculated by simulating the gas exchange process with AVL BOOST software. The same distribution diagram as in the real engine has been considered and the intake pressure has been set to obtain the same average mass flow rate as the measured mass flow rate at each test point. The results are shown in Table 2. 4.3. Turbulence intensity Turbulence intensity (uʹ) is calculated using the simplified K-k model [22] applied to a closed system. This model has as initial conditions at intake valve closing time an average kinetic energy of the mixture KIVC and a turbulent kinetic energy kICV dependent on uʹ. Energies evolve according to the equations of the model Eq. (3). dK dt = − P+K ˙ ρ u ρ u dk dt =P−m ε +k ˙ ρ u ρ u k=3 2muʹ2P=0.3307ct 3 2 √K clV Ap uʹ ε =uʹ3 clV Ap (3) Where Pis the turbulent kinetic energy production term from the average kinetic energy of the mixture, ε is the turbulent kinetic energy dissipation term, Vis the combustion chamber volume and Apis the piston area. Initial kinetic energies KIVC and kICV have been chosen to match those obtained with the gas exchange model developed in AVL BOOST at the instant of intake valve closing. The model has two tuning parameters, ct=0.6 which is a constant affecting turbulence production and cl=0.3 affecting turbulence dissipation. Both parameters have been tuned following the AVL BOOST guidelines so that uʹat top dead centre has values between one and two times the mean linear piston speed. 4.4. Flame properties 4.4.1. Laminar flame front thickness The laminar flame front thickness (δl) is defined using the diffusion scaling approach. δl= α u/Slas in [23,24]. The thermal diffusivity of the mixture ( α u) is defined as α u=λT/( ρ ucp)[25]. The procedure proposed by Mathur et al. [26] is followed to calculate the thermal conductivity of the mixture λT. The cpof the mixture has been calculated following a molar fraction weighted formulation by calculating the cpof each species with the correlations of NIST [27]. 4.4.2. Effect of instabilities on flame front thickness Instabilities cause an increase in flame front thickness which is quantified by the parameter kδlintroduced in [14] Eq. (4). kδl=exp(2 πω DL −L σ (1+ ω DL)( σ + ω DL) σ + ( σ +1) ω DL 4 π 2 LI)(4) This term includes hydrodynamic and thermo-diffusive instabilities. The term ω DL is only dependent on the density ratio σ . Furthermore, Lis the Markstein number, which depends on the effective Lewis, Zel’dovich number and σ [28] Eq. (5). L=δl[ σ ln σ σ −1+β(Leeff −1) 2( σ −1)∫ σ 1 lnx x−1dx ](5) 4.4.3. Zel’dovich number Zel’dovich number is defined as β= (Ea(T0 b−Tu))/(RT02 b). Where T0 b is the adiabatic flame temperature and Eais the activation energy. The activation energy can be calculated according to Eq. (6) [29]. Fig. 2. 1/4 cylinder head views with one of the section planes used to calculate the intersection of the burned volume with the cylinder head. Fig. 3. Autodesk Inventor representation showing the interaction of the flame front with the combustion chamber walls. P. Gabana et al. Fuel 381 (2025) 133563 5 Ea= − 2R⎡ ⎢ ⎢ ⎢ ⎣ ∂ ln( ρ uSl) ∂ (1 T0 b)⎤ ⎥ ⎥ ⎥ ⎦p (6) 4.4.4. Effective Lewis Effective Lewis is defined as Leeff = α u/Di. The mass diffusivity (Di) is a property of a species in the mixture calculated with the Blanc’s law presented in [25]. In this work binary mass diffusivities have been calculated according to [25]. However, when there are mixtures of various fuels, there are different methodologies to calculate Leeff . In [23] various methodologies for the calculation of the effective Lewis of several hydrocarbon and hydrogen mixtures are analysed. The effective Lewis have been calculated according to the heat release-based approach Eq. (7) using the results obtained in [23] and the hydrogen concentrations in the experimental data set of this work. Leeff =1+q1(Le1−1) + q2(Le2−1) q1+q2 (7) Where q1and q2correspond to the energy released by each fuel in the mixture and Le1and Le2are the Lewis numbers of each fuel species. 5. Results 5.1. Thermodynamic diagnostic model results From the pressure signal, mass fraction burned (MFB) is calculated. From MFB, CA timings can be calculated. Fig. 4 (a) shows the average CA timing from the start of combustion to MFB =0.1 and from MFB =0.1 to MFB =0.9. CA timings are calculated from mass fraction burned (MFB) and pressure signal. CA0-90, CA0-10 and CA10-90 are considered to quantify the duration of the total combustion process, the first phase of combustion and the main or second phase of combustion, respectively. The addition of hydrogen clearly decreases the duration of the combustion process for the same air to fuel ratio (Fig. 4 (a)). In G, the first phase of combustion has a longer duration than the second phase. In M the opposite occurs (CA10-90 <CA0-10) for all λvalues lower than 1.9. For both cases (G and M), when λincreases, the share of the first phase of combustion with respect to the total combustion time increases, in agreement with the results of [30,31]. When comparing the results of G with Mfor λvalues between 1 and 1.2, it is concluded that the addition of hydrogen to the fuel leads to a reduction in the first phase of combustion (CA0-10) of between 40 % (λ=1) and 50 % (λ=1.2). However, for the second phase of combustion the reduction ranged from 8 % (λ= 1) to 16 % (λ=1.2). Therefore, the addition of hydrogen overall reduced the combustion duration but affecting more to the first phase than the second phase of combustion. CA0-10 is divided into CA0-1 and CA1-10, Fig. 4 (b). CA0-1 includes the growth of the flame kernel and the transition to turbulent regime [32]. CA0-1 is larger than CA1-10 for all experiments. Both CA0-1 and CA1-10 increased with λfor G and M fuels, but CA0-1 increased much more than CA1-10. Therefore, most of the increase in CA0-10 with increasing λis due to the increase in CA0-1. Fig. 5 (a) shows Sdvalues versus Rffor 200 cycles of the G1.1 experiments. Comparing the Sdresults for each cycle for the same flame radius enables to analyse the cycle-to-cycle fluctuations and to determine the averaged combustion speed when the flame front passes through a certain position. Fig. 5 (b) shows the coefficient of variation (CoV) as a function of flame front radius. The dispersion for all radii between 10 and 40 mm is between 9 % and 14 %. It should be noted that all the experiments have been carried out with the engine running at the same engine speed and, therefore, the turbulence levels are very similar for all experiments. Fig. 6 (a) shows Rfas a function of MFB averaged for all experimental conditions. The mean MFB value is also plotted and highlighted, obtaining very similar values for all experimental conditions. Fig. 6 (b) shows the cycle averaged combustion speed (Sd) at the same flame front position (same Rf) for 200 cycles at the 14 experimental conditions. On the bottom horizontal axis an approximate value of the mass fraction burned (MFB) is presented. The values for flame front radii Rf<10 mm correspond to MFB <1%. The combustion speed values in this zone (Fig. 6 (b)) are not reliable because there is not enough resolution in the pressure signal derivative. Thus, it is not possible to obtain reliable information about the kernel growth. Therefore, the results of Fig. 6 (b) are analysed in more detail for flame front positions larger than 10 mm where MFB >1 %. In the Rfrange from 10 mm to 20 mm (CA1-10 approximately) the Sd of Gare lower than those of M. A large influence of λon Sdfor λvalues between 1 and 1.2 is not obtained in either of the two fuels. In Mwhen the value of λexceeds 1.2, in the Rfrange between 15 and 20 mm, the influence of λis slightly appreciated, so that Sddecreases when λincreases. These trends are shown in CA1-10 (Fig. 4 (b)). CA1-10 values do not significantly vary with λ, being Gvalues larger than those obtained for M. Increasing λin Mexperiments slightly increases CA1-10. This combustion phase takes place in an angular interval between 6 and 9 crank angle degrees. In the second phase of combustion quantified by CA10-90 (approximately Rf>20 mm and Rf<40 mm), the influence of λon Sdfor both Fig. 4. (a) Duration of the first phase of combustion (CA0-10), the second phase of combustion (CA10-90) and the total combustion duration (CA0-90) as a function of λfor Mand Gexperiments. (b) Values of CA0-1 and CA1-10 as a function of λfor Mand Gexperiments. P. Gabana et al. Fuel 381 (2025) 133563 6 Mand Gis observed. In this second phase, combustion speeds of G cease to be lower than those of M.Sdtrends are reflected in the CA10-90 values in Fig. 4 (a), the higher the Sdthe lower the CA10-90 values. As combustion evolves, Sdvalues for Ggrow faster than those of M results. The point at which the flame front reaches the piston, approximately Rf=14 mm, is also observed. This point is very similar for all the experiments, since the combustion process takes place at crank angle times when the piston is close to top dead centre. Fig. 6 (c) shows Slcycle averaged at each experiment as a function of Rf. Each Slhas been calculated with Eq. (2) as a function of the thermodynamic conditions of the mixture at each instant. Slvalues increase with decreasing λ. The addition of hydrogen makes Slvalues with M larger than with Gfor the same λvalues. The difference that exists between Slvalues are weakly reflected in Sdin the first phase of combustion Rf<20mm. In the second phase, Slvalues are ordered with λ, the higher λthe lower Sl. At Rf=30mm an increase in Slaround 550 % is reflected in an increase around 30 % in Sd. Turbulence influence increases strongly with decreasing Sl. The Flame Speed Ratio is defined as FSR =Sd/Sl.Fig. 6 (d) shows ln(FSR −1)cycle averaged as a function of Rffor all experiments. It is observed that FSR depends on Rfbut more on λ. The FSR trend is explained by the fact that Slvaries largely with λwhile Sdtakes values in the same order for all λvalues. According to Gülder [33],FSR is primarily dependent on the ratio uʹ/Sl. The uʹvalue is very similar at all experiments since the engine speed is the same. Therefore, the variation of FSR between all experiments is due to the variation of Sl. Sdvalues for different Rfare shown as a function of λ,Fig. 7. For λ values lower than 1.2 and Rf≤20mm in both Gand Mthere is no clear trend with λ. For higher Rfvalues, the trend of Sdwhen increasing λis as expected from the trend of Sl. In the case of Mwith λ>1.2 the trends are decreasing with λat all Rfvalues. The trend of decreasing Sdwith increasing λis larger as Rfincreases. 5.2. Prediction of flame speed ratio (FSR) Gülder [33] proposed a relationship between the turbulent (St) and the laminar combustion speed (Sl) depending on uʹ/Sland Li/δl[28]. This relationship describes the free-field combustion behaviour, Eq. (8). Fig. 5. (a) Sdas a function of Rffor all cycles of G1.1. (b) Coefficients of variation of Sdas a function of Rffor all experimental points. Fig. 6. (a)MFB cycle averaged as a function of Rfand in thicker line the average of all experimental conditions. (b) Sdcycle averaged as a function of Rf. (c) Slcycle averaged as a function of Rf. (d) FSR for all experiments cycle averaged as a function of Rf. P. Gabana et al. Fuel 381 (2025) 133563 7 St Sl −1=k0(uʹ Sl)a(Li δl)b (8) For the case of internal combustion engines, Gim´ enez et al. [14] proposed a new correlation which considers the effect of instabilities increasing the apparent thickness of the flame front (δl) through a laminar flame front thickness multiplier kδl. Gim´ enez et al. [14] also added a new term to quantify that the flow is confined by the walls of the ICE combustion chamber (confined combustion). It was suggested that having a confined combustion produces an additional turbulence that is taken into account with a new term called as the entrainment velocity (Se), defined as the speed of the reactants close to the flame front. This speed, close to the walls, produces shear stresses that cause further turbulence. The amount of turbulence that is generated and affects the combustion speed depends on the length of the intersection between the flame front and the combustion chamber walls. In this work it is referred to as the flame front perimeter Pf.Considering Sethe intensity of the new turbulence generated, the new correlation term includes the relation between the expansion speed Seand the laminar combustion speed Sl and also the influence of the perimeter of the flame front nondimensionalised with the integral scale LiEq. (9). St Sl −1=a(uʹ Sl)b(Li δlkδl)c(SePf SlLi)d (9) According to Gim´ enez et al. [14] if mean properties over the whole combustion chamber are considered, Secan be calculated with Eq. (10). Se=1 Af[dV dt (Apu Ap −Vu V)+Vu ˙ mbb ]+St((γ−1)muHp γpV )+(γ−1) γpAf(˙ Qw Vu V− ˙ Qu)(10) The correlation proposed by Gim´ enez et al. [14] has the drawback that Sedepends on Stso it is not an explicit expression. In this work, an iterative new method to calculate Stis used. The method consists of making a hypothesis of Stand substituting it in Eq. (10) to calculate Se and then the new Stwith Eq. (9). The process is repeated iteratively until both Stvalues match. In [14] the coefficients were fitted using multilinear regressions. However, the current methodology does not allow to fit the coefficient aand the exponents b,cand dto the experimental data using the method proposed in [14]. In this work, the coefficients have been fitted using a new optimisation method. The method is based on using a genetic algorithm integrated in Matlab in which the function to be optimised is the mean square error between Sdand the predicted combustion speed St. For 12 specific values of the flame front radius between 10 and 37 mm, the average values of all variables involved in the correlation have been calculated over 200 cycles, resulting in 12 different values of all variables at the 14 experimental conditions. The results from the fit are shown in Fig. 8. Table 3 shows the correlation coefficients obtained in [14], the coefficients obtained in this work and the coefficients of the fit of the data from [14] using the methodology of this work. The exponents of the correlation obtained following the new methodology proposed in this work, where Sedepends on St, are very similar for the two data sets (cases 2 and 3). There is a remarkable difference, especially in the exponent of the third term, d, when it is calculated with the experimental values of Sd(cases 1 and 3). The third term is calculated with Stand it is the term modified in the iterative process. The estimation of uʹgreatly affects the coefficient aof Eq. (9). The value of aobtained in case 2 is lower than in case 3. In this work the values of uʹat top dead centre divided by the piston mean linear velocity are larger than in [14]. An underestimation of u’in [14] or an overestimation in this work would justify the differences in the value of a. Even with very different combustion chamber geometries (pent-roof cylinder head vs flat cylinder head) the correlation exponents take very similar values in cases 2 and 3. Fig. 9 shows the Borghi-Peters diagram [34,35] and the location of the experimental conditions for four flame front positions (15, 22, 30 and 37 mm), under two laminar flame front thickness hypothesis δland δlkδl. Three iso-FSR lines have been added (FSR =2, FSR =10 and FSR = 30) into Fig. 9. These iso-FSR lines have been calculated by using the expression of Gülder [33] corresponding to free field combustion. Considering δlas the flame front thickness, when λincreases, as Sl decreases and uʹhas similar values, the points move from the zone with Lk>δlto the zone Lk<δlby changing the combustion regime. Therefore, if it is assumed that the flame front thickness is increased by the instability (flame front thickness as δlkδl) the combustion regime in all experimental conditions is in the thin reaction zone instead of close to Fig. 7. Sdas a function of λfor different Rf. Fig. 8. Experimental data of ln(Sd /Sl −1) versus values obtained from the correlation. The plot includes12 flame front positions for Mand G. P. Gabana et al. Fuel 381 (2025) 133563 8 the corrugated flames zone. It is observed that the test points would be located approximately between FSR =10 and FSR =30 [33] when the flame front thickness is δl. If the influence of instabilities in the flame front thickness (δlkδl) is taken into account, the points are shifted to the left with lower FSR values. 5.3. Understanding the effect of the correlation terms on FSR prediction First and third terms of Eq. (9) quantify the relationship between the turbulence intensity in the combustion chamber and the laminar combustion speed. The turbulence intensity can be caused by the fluid velocity in the combustion chamber or by the expansion of the reacting gases. The second term reflects the influence of the ratio of the characteristic turbulence dimension Lito the flame front thickness modified by instabilities. Fig. 10 shows the three combustion terms as a function of the flame front radius for all experimental conditions. The trend of uʹ/Slwith Rfis shown in figure Fig. 10 (a). The trend is always decreasing as Rfincreases since the turbulence intensity decreases with time due to its dissipation and Slincreases due to temperature increase. In Fig. 10 (c) the trend of the term (SePf)/(SlLi)with Rfis shown. This term is the one that provides the shape of Stas a function of Rf. As mentioned above, this term is attributed to the turbulence generated by the unburned gas speed with respect to the combustion chamber walls. The trend with respect to λis the same as uʹ/Sl, mainly due to the decrease of Sl. In first and third terms of Eq. (9), the representative curves of each λ do not intersect when changing Rf, so they can hardly explain the different trend when varying the radius of Gand MFig. 4 (a). Fig. 10 (b) shows that Li/(δlkδl)term does not maintain the same trend against Rf for all λvalues. The value of this term is not ordered as a function of λfor a fixed Rfas in the other two correlation terms. For the same Rfvalue, the highest values of this term are reached at λvalues of about 1.5 in M. The term takes similar values for λ=1 and for λ=2 in Mfor Rfvalues lower than 20 mm. Subsequently, at values of Rf>20mm the term takes larger values in M1 than in M2. Gvalues have trends with Rfsimilar to those in M1, but with smaller values. In Mresults, the higher the λ, the higher the value of the term. The evolution of this term with Rfcan explain the different behaviour of Gresults compared to M. Fig. 11 shows predicted FSR as a function of Rffor 4 test points, the experimental FSR are also plotted. The predicted trends with fuel type, λ and with Rfare acceptable especially in the central zone (Rf>10mm and Rf<37mm). Nothing can be stated about the start of combustion since the experimental results cannot be calculated. Regarding the end of combustion, the prediction does not decrease at the same rate as the experimental data. In this zone the effects of quenching can be proportionally important with respect to the rate of fuel being burned. In predictive models where this correlation is used, it is necessary to complete them with kernel growth and flame quenching models. Table 3 Fitting results with different data and methodology. Case Data Methodology a b c d 1[14] [14] 0.114 0.86 0.29 0.4 2 This work This work 0.124 0.8 0.21 0.32 3[14] This work 0.21 0.85 0.26 0.31 Fig. 9. Borghi-Peters diagram and experimental results for 4 flame front positions. The efect of considering the flame front thickness as δl(points on the right) and as δlkδl(points on the left) is shown. Fig. 10. Cycle averaged values of the correlation terms Eq. (9) as a function of Rffor all experimental conditions. (a) uʹ/Sl. (b) Li/(δlkδl). (c) (SePf)/(SlLi). P. Gabana et al. Fuel 381 (2025) 133563 9