**2. Engine System**

The engine used for the experimental activity is a twin-cylinder turbocharged spark ignition engine, equipped with two port injectors, one for each cylinder, to supply the gasoline just upstream of the intake valves. Engine's main characteristics are reported in the following Table 1. It is provided by a conventional pent-roof combustion chamber, a centered spark plug and a standard ignition system. Each cylinder presents 4 valves, two intake valves and two exhaust valves. An electro-hydraulic Variable Valve Actuation (VVA) module is mounted on the intake side, allowing for a flexible control of the lift profile. This device provides the actuation of both Early Intake Valve Closure (EIVC) and Full Lift valve strategies. On the exhaust side, a fixed valve lift strategy is employed. Engine boosting is realized by a small waste-gated turbocharger. The compressor operating domain bounds the engine performance because of surge and choke phenomena and of a maximum allowable rotational speed (255,000 rpm). In particular, at low speeds, the boost pressure is limited to avoid the compressor surging occurrence [17]. An additional constraint is imposed by the manufacturer to the maximum boost level, achieved at high speeds, to guarantee the mechanical integrity of the intake plenum.

**Table 1.** Engine characteristics.


The engine is designed following the so-called "downsizing" concept. At high loads the engine works under knock-limited conditions and this requires to delay the combustion phasing (50% of mass fraction burned—MFB 50%), especially at lower speeds where a higher knock tendency usually occurs. On the other hand, at high speeds the A/F mixture has to be particularly enriched to maintain the Temperature at Turbine Inlet below a certain maximum allowable level. These control strategies, although mandatory for the examined engine, greatly penalize the fuel consumption at high loads. For this reason, the base engine configuration is properly modified for research purpose through the installation of an external low-pressure (LP) exhaust gas recirculation (EGR) circuit as depicted in Figure 1.

**Figure 1.** Engine layout with LP EGR circuit.

EGR rate and temperature are controlled through the throttle EGR valve and the cooler device, respectively. An enhanced cooling system, requiring a water cooled heat exchanger, is utilized to cool down the recycled gas. The LP EGR system was preferred over a high-pressure (HP) device, because it allows to avoid the pressure fluctuations and the counter-flow within the EGR circuit, due to the turbocharger dumping effect.

As known, EGR is essentially considered for knock control purposes at high loads, although certain EGR-related efficiency advantages are also obtained at part loads, due to the engine de-throttling.

The engine at test bench is also equipped along the exhaust line with a Three-way Catalytic converter (TWC) to guarantee the pollutant emissions (HC, CO and NOx) abatement, provided that the air/fuel ratio window very close to the stoichiometric value is maintained.

The exhaust system was opportunely modified to extract a portion of the exhaust gas from a single cylinder. In particular, a thin metallic tube is inserted inside of the exhaust manifold, through a properly realized hole, and oriented towards the exhaust side of one cylinder. The internal end of the above tube is located very close the exhaust valves. In this way, the exhaust gases of the selected cylinder (Cyl #1 in the following) are extracted and externally derived to measure the individual cylinder-out noxious emissions.

#### **3. Experimental Setup and Test Procedure**

As aforementioned, the tested engine is provided by an electro-hydraulic variable valve actuation (VVA) device. In spite of the VVA potentials, in the experimental tests the load is adjusted only acting on the throttle valve opening or waste-gate valve opening, without modifying the intake valve lift profile which is set to the "Full Lift" configuration. The intake air is constantly supplied to the engine at 293 ± 1 K by an air conditioning unit. Each engine cylinder is equipped with piezo-quartz pressure transducer (accuracy of ±0.1%) to detect the in-cylinder pressure signal.

The instantaneous pressure signals are acquired by the AVL INDICOM over 270 consecutive pressure cycles for the combustion analysis, assuming a polytropic thermodynamic process. A resolution of 0.1 CAD within the angular window between −90 and 90 CAD AFTDC is chosen, while outside this angular interval the sampling resolution is set at 1 CAD.

An automatic post-processing tool, based on a thermodynamic model, computes the ensemble averages of in-cylinder pressure and burn rate profiles and the combustion characteristics data. Furthermore, the boost pressure and the upstream turbine pressure are acquired through piezo-resistive low pressure indicating sensors located at the compressor outlet and at the outlet of the exhaust manifold, respectively. The engine is also equipped with thermocouples to monitor intake and exhaust temperatures, with particular attention to the control of the turbine inlet temperature in order to avoid unacceptable levels for the turbine blades.

A prototype driver, managed within LabView environment, is capable to switch from the commercial ECU to an external user control of the main engine variables: fuel injection timing and duration, boost pressure, spark timing, external EGR flow and temperature. EGR-related variables are managed acting on valve opening and cooler device within the EGR circuit.

Exhaust gas emissions are sampled upstream of the three-way catalyst. An ultra-violet gas analyzer (ABB UV Limas 11) measures NOx. A cold extractive IR gas analyzer (ABB URAS) detects CO, CO2 and O2 while a FID analyzer (Siemens, Milano, Italy) is used for THC. A non-dispersive infrared sensor measures the CO2 concentration at the inlet, downstream of the throttle valve to monitor the EGR rate. The sensor error is set at 0.07% for the measured CO2 concentration.

The engine is tested in steady-state conditions at two different speeds (1800 and 3000 rpm) and various loads (from low to medium/high levels), also including increasing external EGR rates. In particular, EGR is modified by keeping constant the engine IMEP and the relative A/F ratio, λ. This means that at increasing the EGR rate the IMEP is restored to the nominal value with engine boosting; in a similar way, λ is maintained to the stoichiometric value by a modulation of gasoline injection duration. These operating points are selected because they are representative of the WLTC cycle for a segmen<sup>t</sup> A vehicle equipped with the engine under investigation. Performance in these points, such as fuel consumption and CO2 emissions, are relevant for the vehicle homologation. Various operating parameters are collected during the experiments, including torque/power, air flow rate, fuel flow rate, boost and exhaust pressures and temperatures, exhaust emissions, etc.

The overall considered operating conditions are reported in Table 2. Test grid collects 31 points which were gathered into 7 groups, each one characterized by the same nominal load level and rotational speed. For a single group, a label is defined (Table 2), referring to the nominal net IMEP and speed, which is utilized in the following figures. The table also shows the measured values of EGR rate, relative A/F ratio, λ, and the spark advance (SA). This last is selected to realize the Maximum Brake Torque (MBT) condition at knock free operations, while at high/medium loads the SA is chosen to operate at knock borderline (knock limited spark advance, KLSA). In particular, knocking is monitored by a knock sensor, installed between the two cylinders on the block, which automatically cuts out the power output and avoids the engine operation under severe knock conditions.

**Table 2.** Overall set of experimental engine points.


The measured net IMEP (Table 2) is properly derived from the acquired in-cylinder pressure traces over the entire engine cycle and represents the difference between the gross IMEP (in-cylinder pressure over compression and expansion strokes) and the pumping mean effective pressure (PMEP), evaluated over the intake and exhaust strokes.

Further parameters are monitored to preserve the safety of the turbocharger and of the entire engine. In particular, the boost pressure is controlled at high loads, acting on the waste-gate (WG) valve opening, to provide proper operation for the port injectors and to avoid the mechanical failures for the intake manifold. The average maximum in-cylinder pressure is constrained to limit the engine mechanical stresses. The mechanical and thermal safety of the turbine also obliges to control the turbocharger speed and the turbine inlet temperature (TIT). Summarizing, the following constraints are taken into account:


## **4. Experimental Results**

As discussed in the introduction section, experiments carried out on the examined engine have highlighted significant differences in the combustion evolution between cylinders, mainly ascribed to a non-uniform effective in-cylinder air/fuel (A/F) ratio [15]. Cylinder-by-cylinder variations are apparent in Figure 2a which shows the experimental incylinder pressure trace and the rate of heat release (ROHR) for two representative operating points (3000@9 and 3000@13, 0% EGR). As aforementioned, the pressure data refer to the ensemble average over 270 consecutive cycles. For both operating conditions, the pressure curves of two cylinders are overlapped in the compression stage, denoting the same air volumetric efficiency, while Cyl #2 provides higher pressure peak and combustion rate than Cyl #1. The same behavior is found at each investigated operating point, suggesting a systematic difference in pressure cycles. Further experimental tests on the adopted fuel injection system have highlighted a different fuel supply to cylinders, resulting from a variation in port injectors fuel rates (Figure 2b). The results plotted in Figure 2b were collected through injection tests realized at am-bent air pressure with fuel injection system removed from the engine. A drive signal to injectors was used to simulate the injection timing; a certain number of injection timings were considered to measure the gasoline mass per stroke delivered by injectors. In particular, consistently with the in-cylinder pressure traces shown in Figure 2a, the injector corresponding to Cyl #2 provides a higher fuel mass flow rate. Interestingly, no difference was found also when the relative position of two injectors was switched, indicating the fuel rail geometry as responsible of the cylinder-by-cylinder variation. Indeed, the gasoline rail mounted on the examined engine at test bench represents a non-optimized prototype geometry. The optimized version of rail, mounted on the commercial SI engine, does not exhibit differences in the gasoline mass flow rates between port injectors. All the previous considerations underline that rail geometry is main responsible for the experimental cylinder-by-cylinder variation analyzed in this research activity.

**Figure 2.** Comparison of experimental in-cylinder pressure traces and burn rates profiles at 3000 rpm, two IMEP levels (9 and 13 bar) and EGR = 0% (**a**). Differences of the injected gasoline mass per stroke between port injectors (**b**).

Figure 3 shows the relative air–fuel ratio (λ) in Cyl #1, Cyl #2 and the engine exhaust at different engine operating points and EGR rates. λ of Cyl #1 was estimated from the composition of the exhaust gases of that cylinder through the carbon balance of species. Similarly, the overall engine λ was verified from the overall exhaust gas composition. This allows to indirectly evaluate the λ value of Cyl #2 under the hypothesis of an equal volumetric efficiency for both cylinders, as suggested by the comparison of pressure traces in the compression stage. Even though the engine conditions are close to stoichiometric, lean (λ between 1.03 and 1.08) and rich (λ between 0.96 and 0.98) mixtures are obtained in Cyl #1 and Cyl #2, respectively. Considering that the engine does not provide an individual set of injection and combustion phasing for each cylinder, the same spark advance and injection parameters are chosen for both cylinders, optimizing the engine load. As a consequence, a difference in the single cylinder load is found at each operating condition.

**Figure 3.** Relative air–fuel ratio (λ) in Cyl #1, Cyl #2 and overall engine at different operating points as a function of EGR rate.

The individual IMEP levels for both cylinders at different operating points are shown in Figure 4 as a function of the EGR rate. As expected, the larger fuel amount introduced in the Cyl #2 results in a higher IMEP at each investigated condition. Furthermore, the adoption of the same spark advance for both cylinders penalizes the combustion phasing of Cyl #1 which results slightly late compared to the MBT. As concerns the trend versus the EGR rate, it is worth noting that the measurements are performed keeping constant the overall engine load; hence, the single cylinder IMEPs are almost constant at varying the exhaust-gas recirculation.

In Figure 5, the EGR-related trend of the engine ISFC is presented for the measured operating points. As expected, a reduced ISFC is observed at low and medium/high IMEPs by increasing the EGR content. At low loads, the ISFC benefits are ascribed to the EGR-induced reduction in the pumping losses promoted by the engine de-throttling; at medium/high loads, the ISFC advantages are obtained thanks to the EGR capability in reducing the knock tendency [18].

**Figure 4.** Indicated load in Cyl #1 and Cyl #2 at different engine operating points as a function of EGR rate.

**Figure 5.** Engine Indicated Specific Fuel Consumption at different operating points as a function of EGR rate.

Figures 6 and 7 show the characteristic angle of MFB 10% and the combustion core duration (MFB 10–50%) for both cylinders at different operating conditions. Looking at the early stage of combustion, at each operating condition, a slight delay in MFB 10% (Figure 6) of Cyl #1 compared to Cyl #2 is evident: the richer mixture in Cyl #2 provides a higher flame speed, advancing the MFB 10%. As shown in Figure 6, this gap increases at higher EGR rates; as is well-known, the charge dilution slows down the combustion rate, enhancing the effect of the equivalence ratio on the flame speed. The cylinder-to-cylinder variation in MFB 10% does not provide significant differences in the core combustion duration as reported in Figure 7. Even though the combustion duration in Cyl #1 results always prolonged compared to Cyl #2, the maximum difference in MFB 10–50% is less than 1.5 CAD at 1800@10 with 14.4% of external EGR. Regarding the effect of charge dilution on the combustion duration, the EGR prolongs the MFB 10–50% at each operating condition.This trend is more evident at higher engine speed, due to the reduction in cycle length, leaving a shorter period available to complete the combustion process. Further increase in combustion duration against the EGR is found at lower engine load, because of the greater impact of internal EGR which contributes to increase the overall in-cylinder residual content at decreasing the load.

**Figure 6.** MFB 10% in Cyl #1 and Cyl #2 at different operating conditions as a function of EGR rate.

**Figure 7.** MFB 10–50% in Cyl #1 and Cyl #2 at different operating conditions as a function of EGR rate.

#### **5. Description of Modeling Approach**

The engine used in this paper is outlined in a 1D model by employing the GT-Power commercial code. The whole engine is represented in sub-components, including cylinders, intake/exhaust pipes and turbocharging system. In particular, the flow inside the intake/exhaust ducts is reproduced by solving the conservation equations of mass (1), energy (2) and momentum (3) reported below:

$$\mathbf{\dot{d}m/dt} = \boldsymbol{\Sigma}\dot{\mathbf{m}},\tag{1}$$

$$\mathbf{d(m\cdot e)/dt} = -\mathbf{p} \cdot \mathbf{d} \mathbf{V}/\mathbf{dt} + \Sigma(\mathbf{rih}) - \mathbf{dQ\_w}/\mathbf{dt},\tag{2}$$

$$\mathbf{dm}/\mathbf{dt} = \left(-\mathbf{dp}\cdot\mathbf{A} + \Sigma(\dot{\mathbf{m}}\cdot\mathbf{u}) - \mathbf{4}\cdot\mathbf{C}\_{l}\cdot\mathbf{p}\cdot\mathbf{u}\cdot\|\mathbf{u}\|\cdot\mathbf{A}\cdot\mathbf{dx}/2\cdot\mathbf{D} - \mathbf{A}\cdot\mathbf{C}\_{p}\cdot\mathbf{p}\cdot\mathbf{u}\cdot\|\mathbf{u}\mid/2\}/d\mathbf{x},\tag{3}$$

where m, m, V, p, e, h and u are mass flow rate, mass, volume, pressure, internal energy ˙ per unit mass, enthalpy per unit mass and velocity at boundary. Q w is the heat exchanged through the walls, A is the flow area, D is the equivalent diameter. Cp and Cf are the pressure and friction losses coefficients, while dx and dp represent the discretization length and the pressure differential across dx. Each cylinder is treated as a 0D volume. In this case, the scalar variables (e.g., pressure, temperature, internal energy, etc.) are assumed to be uniform over the volume. However, during the combustion process the cylinder volume is divided in two zones (burned and unburned regions), where the mass and energy conservation equations are solved. In the closed valve period of single cylinder (i.e., from IVC up to EVO), the energy equation is detailed for burned and unburned zones as reported by relations (4) and (5):

$$\mathbf{d}(\mathbf{m}\_{\mathbf{u}} \cdot \mathbf{e}\_{\mathbf{u}})/\mathbf{dt} = -\mathbf{p} \cdot \mathbf{d} \mathbf{V}\_{\mathbf{u}}/\mathbf{dt} - \mathbf{d} \mathbf{Q}\_{\mathbf{w},\mathbf{u}}/\mathbf{dt} - \mathbf{h}\_{\mathbf{u}} \cdot \mathbf{dm}\_{\mathbf{b}}/\mathbf{dt},\tag{4}$$

$$\mathbf{d}(\mathbf{m}\_{\mathrm{b}} \cdot \mathbf{e}\_{\mathrm{b}})/\mathrm{dt} = -\mathbf{p} \cdot \mathrm{dV}\_{\mathrm{b}}/\mathrm{dt} - \mathrm{d}\mathrm{Q}\_{\mathrm{w}/\mathrm{b}}/\mathrm{dt} + \mathrm{h}\_{\mathrm{u}} \cdot \mathrm{dm}\_{\mathrm{b}}/\mathrm{dt},\tag{5}$$

where dmb/dt represents the burning rate term evaluated through the turbulent combustion sub-model discussed in the following.

The turbocharging operation is simulated by a "map-based" approach. Turbine and compressor maps, supplied by the manufacturer, are employed within the standard turbine and compressor objects, respectively. Flow permeability of the single cylinder head is taken into account through the measured steady flow coefficients, both under forward and reverse flow conditions. Standard injector objects, available in the 1D code, are adopted for gasoline injection. During port injections, it is assumed that 30% of the total liquid mass vaporizes immediately upon the injection event, as advised for a typical liquid injector in GT-Power user manual; details about the spray evolution and wall film formation are not taken into account. In order to enhance the model capability, it is integrated with refined phenomenological in-cylinder sub-models of turbulence, combustion, heat transfer and pollutant emissions. A detailed description of these sub-models is reported in the next subsection.

#### *5.1. In-Cylinder Sub-Models and Tuning: Combustion, Turbulence and Pollutant Emissions*

In-cylinder 0D sub-models were integrated within the 1D engine model to refine the prediction of turbulent combustion and of the exhaust emissions. The combustion process is modeled by the fractal approach [19–21], where the burning rate is written as follows:

$$\mathbf{d}\mathbf{m}\_{\rm b}/\mathbf{dt} = \rho\_{\rm u} \cdot \mathbf{S}\_{\rm L} \cdot \mathbf{A}\_{\rm T} = \rho\_{\rm u} \cdot \mathbf{S}\_{\rm L} \cdot \mathbf{A}\_{\rm L} \cdot (\mathbf{L}\_{\rm max}/\mathbf{L}\_{\rm min})^{(\rm D\_3 - 2)},\tag{6}$$

where ρu is the unburned gas density, AL and AT the area of the laminar and turbulent flame fronts, respectively, and SL the laminar flame speed (LFS). Lmin and Lmax are the minimum and maximum flame wrinkling scales and D3 is the fractal dimension. Lmin represents the Kolmogorov scale evaluated under the hypothesis of isotropic turbulence, while Lmax is proportional to a characteristic dimension of the flame front, i.e., the flame radius in the present study. D3 is computed by the equation proposed in Reference [22], as a function of the ratio between the in-cylinder turbulence intensity and LFS. The latter is evaluated by using an "in-house" developed correlation, based on 1D LFS computations via chemical kinetics solver (CANTERA), accounting for the in-cylinder thermodynamic state (pressure p and temperature T), equivalence ratio, Φ, and EGR dilution. In particular, the adopted correlation reported in Equation (7) presents the well-known expression of the power law:

$$\mathbf{S}\_{\rm L} = \mathbf{S}\_{\rm L0} \cdot (\mathbf{T}/\mathbf{T}\_{\rm ref})^{\alpha} \cdot (\mathbf{p}/\mathbf{p}\_{\rm ref})^{\beta} \cdot \mathbf{EGR}\_{\rm factor} \tag{7}$$

In Equation (7), SL0 is the flame speed at reference conditionsT=Tref and p = pref, depending on the fuel sensitivityS=RON-MON and Φ; EGRfactor is a reduction term for LFS that accounts for the presence of residual gas in the unburned mixture. Exponent α depends on Φ, T and p, while exponent β includes the dependencies on Φ and S. It is the case to underline that the above LFS correlation was developed considering a maximum EGR percent mass of 20%; therefore, the employed LFS correlation shows the capability to cover the range of examined EGR rates during the experiments (Table 2). The laminar flame area, AL, is computed by an automatic procedure implemented in a CAD software and processing the actual 3D geometry of the combustion chamber. The estimation of Lmin, and D3 is based on the K–k–T turbulence sub-model, extensively reported in Reference [23]. The governing equations for the mean flow kinetic energy, K, the turbulent kinetic energy, k, and the tumble momentum, Tm, are as follows:

$$\mathbf{d}(\mathbf{m} \cdot \mathbf{K})/\mathbf{dt} = (\dot{\mathbf{m}} \cdot \mathbf{K})\_{\mathrm{inc}} - (\dot{\mathbf{m}} \cdot \mathbf{K})\_{\mathrm{out}} - \mathbf{f}\_{\mathrm{d}} \cdot \mathbf{m} \cdot \mathbf{K}/\mathbf{t}\_{\mathrm{Turn}} + \mathbf{m} \cdot \mathbf{K} \cdot \dot{\boldsymbol{\Psi}} / \boldsymbol{\rho} - \mathbf{P},\tag{8}$$

$$\mathbf{d}(\mathbf{m}\mathbf{k})/\mathbf{d}\mathbf{t} = (\mathbf{\dot{m}}\cdot\mathbf{k})\_{\mathbf{\dot{m}}\mathbf{c}} - (\mathbf{\dot{m}}\cdot\mathbf{k})\_{\mathbf{out}} + (2\cdot\dot{\Psi}/3\cdot\mathbf{\rho})\cdot(\mathbf{m}\cdot\mathbf{k} - \mathbf{m}\cdot\mathbf{v}\_{\mathbf{l}}\cdot\dot{\Psi}/\rho) + \mathbf{P} - \mathbf{m}\cdot\boldsymbol{\varepsilon},\tag{9}$$

$$\mathbf{d}(\mathbf{m} \cdot \mathbf{T}\_{\mathrm{m}})/\mathrm{d}\mathbf{t} = (\dot{\mathbf{m}} \cdot \mathbf{T}\_{\mathrm{m}})\_{\mathrm{inv}} - (\dot{\mathbf{m}} \cdot \mathbf{T}\_{\mathrm{m}})\_{\mathrm{out}} - \mathbf{f}\_{\mathrm{d}} \cdot \mathbf{m} \cdot \mathbf{T}\_{\mathrm{m}}/\mathbf{t}\_{\mathrm{Turn}}.\tag{10}$$

In Equations (8)–(10), the first and second terms describe the incoming and the outcoming convective flows through the valves, respectively. The third term for the K and Tm equations is the decay contribution due to the shear stresses with the combustion chamber walls by means of the decay function, fd, and the characteristic time scale, tTum. In Equation (9), third term represents an additive compressibility term, proportional to ( Ϻ /ρ); the latter is also included in the K equation as the 0fourth term. P describes the energy cascade mechanism, which causes the transfer of kinetic energy from the mean flow to the turbulent flow. Finally, the last term in Equation (9) takes into account the viscous dissipation rate of k into heat. The described turbulence sub-model proved to be able to adequately reproduce the in-cylinder turbulence evolution along the whole engine cycle, also sensing the variations in the engine operating conditions, the valve strategies and the cylinder geometrical parameters [23]. The heat-transfer model both for cylinders and exhaust subsystem is considered, applying a wall temperature solver based on a finite element (FE) approach. For the in-cylinder heat transfer (gas-to-wall), the Hohenberg correlation is implemented into the 1D code while convective, conductive and radiative heat transfer modes are considered for the exhaust pipes. The Hohenberg correlation here utilized is illustrated in Equation (11):

$$\mathbf{H} = \mathbf{A} \cdot \mathbf{V}^{-0.06} \cdot \mathbf{p}^{0.8} \cdot \mathbf{T}^{-0.4} \cdot \left(\mathbf{v\_{pm}} + \mathbf{B}\right)^{0.8},\tag{11}$$

where V is the volume, p is the pressure, T is the mean gas temperature and vpm is the mean piston speed. The calibration constants, A and B, are calculated by Hohenberg and here used as 130 and 1.4, respectively. The convective heat transfer coefficients for engine coolants (oil and water) are evaluated by simulations of coolant circuits. They are subsequently imposed as an input in the code, assuming a dependency on the engine speed. Cooling boundary conditions (temperatures for oil and water) are introduced in the model according to the levels recommended by the engine manufacturer. Concerning the model tuning, an integrated 3D/1D approach is adopted for turbulence sub-model [23]; model constants are identified to better reproduce the turbulence outcomes of in-cylinder 3D CFD simulations, under motored conditions and for various speeds and valve strategies [23]. Referring to the combustion model tuning, this sub-model includes three tuning constants differently affecting the distinct phases of combustion process. A single set of tuning constants was identified for the considered operating points, and tuning is kept fixed regardless the engine operating conditions. A particular attention is devoted to the modeling of the pollutant emissions. To this aim, "in-house developed" sub-models for the estimation of the main regulated cylinder-out emissions, namely HC, CO and thermal NO, are implemented into the 1D code. The propagation of the cylinder-out noxious species

within the exhaust system up to the inlet section of TWC is also considered in the 1D model. In particular, CO concentrations are computed by solving, in the burned gas zone, a simplified chemical kinetic sub-model based on two-step reactions as reported in Reference [24]. The CO mechanism comprises the following high-temperature CO oxidation reactions (12) and (13):

$$\text{CO} + \text{OH} \leftrightarrow \text{CO}\_2 + \text{H}\_2\tag{12}$$

$$\text{CO} + \text{O}\_2 \leftrightarrow \text{CO}\_2 + \text{O}\_2 \tag{13}$$

For thermal NO emission, a refined approach is here adopted, where the semi-detailed chemical kinetics scheme proposed by Andrae [25] is applied to the burned gas zone to compute the NO concentrations in each cylinder. The kinetic mechanism includes 5 elements, 185 species and 937 reactions. This methodology, although computationally more expensive than the classic Zeldovich scheme, seems to be better in predicting the NO emissions of the single cylinder at varying the thermodynamic conditions, including pressure, temperature, inert content and mixture quality. HC production is simulated by the crevices model and a wall-flame quenching correlation, combined with a thermal-based HC oxidization equation during the expansion stroke. This means that the mixing-based oxidation contribution is neglected by assuming a perfect and instantaneous mixing of unburned charge from crevices with the hot burned gases. In the filling/emptying crevices model, HC species are assumed to accumulate/be released during the in-cylinder pressure rise/decrease phase in/from an arbitrary assigned constant volume, Vcrev, as a fraction of the combustion chamber volume at TDC [26]:

$$
\langle \chi\_{\text{CTV}\_{\text{V}}} \rangle\_{\otimes} = 100 \cdot \text{V}\_{\text{CTV}} / \text{V}\_{\text{cyl}\text{-TDC}\_{\text{V}}} \tag{14}
$$

The volume Vcrev schematizes the crevices in the combustion chamber where the flame front extinguishes and it is controlled in the model by the input parameter, xcrev. The temperature in crevices volume is considered to be the same as the cylinder wall. A constant flow coefficient is attributed to the orifice of the crevices volume both for incoming and outgoing mass flows. Therefore, the pressure evolution within the crevice volume is computed by combining the related mass and gas temperature profiles. For the flame wall-quenching contribution, a simple correlation is here assumed which furnishes an estimation of the total HC arising from the extinguished flame at cylinder walls as a function of the flame quenching distance. This last parameter is assumed inversely proportional to the laminar flame speed, since the flame reaches the wall in laminar conditions. Flame quenching distance varies during the combustion evolution, and, in this study, it is taken towards the end of combustion process (i.e., 90% of mass fraction burned). Indeed, it is quite reasonable to assume that the HC generated by the flame-quenching during the initial combustion phase rapidly diffuse in the burned zone and oxidize. The correlation employed in this study for the quenching-related unburned HC resembles the simple model proposed in Reference [27], and it is shown in Equation (15):

$$\text{HC}\_{\text{quench}} = \text{H}\_1 \cdot \text{\S}\_{\text{L}} + \text{H}\_{2\text{\textdegree}} \tag{15}$$

where δL is the laminar flame thickness, while H1 and H2 are sub-model tuning constants. Finally, for the oxidation of the overall HC level during the expansion stroke the one-step kinetic equation proposed in Reference [28] is adopted. The rate of post-flame HC oxidation is computed according to Equation (16):

$$\text{[dBIC]}/\text{dt} = -6.7 \times 10^{15} \times \text{e}^{-(18735/\text{R} \cdot \text{T} \cdot \text{Inc})} \cdot \text{[HC]} \cdot \text{[O}\_2\text{]} \text{ (p/R} \cdot \text{T}\_{\text{hc}})^2,\tag{16}$$

This rate depends on the HC and molecular oxygen concentrations, the density term (p/RThc) and on the oxidation temperature, Thc. The latter is computed as a weighted average between the in-cylinder mean gas temperature and the wall temperature.
