Numerical Investigation of Dual Fuel Combustion on a Compression Ignition Engine Fueled with Hydrogen/Natural Gas Blends

The present work aims to assess the influence of the composition of blends of hydrogen (H2) and Natural Gas (NG) on Dual Fuel (DF) combustion characteristics, including gaseous emissions. The 3D-CFD study is carried out by means of a customized version of the KIVA-3V code. An automotive 2.8 L, 4-cylinder turbocharged diesel engine was previously modified in order to operate in DF NG–diesel mode, and tested at the dynamometer bench. After validation against experimental results, the numerical model is applied to perform a set of combustion simulations at 3000 rpm–BMEP = 8 bar, in DF H2/NG-diesel mode. Different H2–NG blends are considered: as the H2 mole fraction varies from 0 vol% to 50 vol%, the fuel energy within the premixed charge is kept constant. The influence of the diesel Start Of Injection (SOI) is also investigated. Simulation results demonstrate that H2 enrichment accelerates the combustion process and promotes its completion, strongly decreasing UHC and CO emissions. Evidently, CO2 specific emissions are also reduced (up to about 20%, at 50 vol% of H2). The main drawbacks of the faster combustion include an increase of in-cylinder peak pressure and pressure rate rise, and of NOx emissions. However, the study demonstrates that the optimization of diesel SOI can eliminate all aforementioned shortcomings.


Introduction
Greenhouse gas (GHG) emissions caused by human activities are responsible, thus far, for an average global warming of about 1.0 • C, compared to the pre-industrial era.
In particular, Figure 1 shows an average increase of the global temperature relative to the period 1961-1990 of approximately 0.7 • C [1]. Recently, scientists warned that global warming could rise even up to 4 • C by the end of the century if a carbon-neutral economy is not be achieved by 2050 [2]. As a result, the 2021 edition of the United Nations Annual Climate Change Conference (COP26) recommended a halving of emissions from all economic sectors over the next decade and reaching net zero carbon emissions by the middle of the century [3].
As far as the road transport sector is concerned, many different technologies can be considered for the abatement of GHG emissions at the tailpipe of Internal Combustion Engines (ICEs). Among them, the use of low-carbon fuels and carbon-neutral fuels (biofuels and e-fuels), present a very promising solution to strongly reduce carbon dioxide (CO 2 ) emissions, while maintaining all practical advantages of conventional powertrains, in comparison to Battery Electric Vehicles (BEVs). Moreover, running ICEs on hydrogen (H 2 ), blended with other fuels or alone (see references [4,5], respectively) yields a further advantage: since H 2 does not contain carbon, its combustion does not produce carbon monoxide (CO), CO 2 , Unburnt Hydrocarbons (UHC) and soot. The only critical pollutants The fundamental role of H 2 in the future mobility scenario, as well as in many other economic sectors, is highlighted, as an example, by the "EU Hydrogen Strategy" [6], included in the framework of the "European Green Deal" [7], and reaffirmed in the "NextGenerationEU" recovery plan [8]. The European Union clearly intends to promote environmentally-friendly technologies focused on the use of green hydrogen.
Therefore, in addition to fuel cells and BEVs, ICEs can continue to represent a fundamental technology for vehicle propulsion, provided that their main drawbacks (i.e., CO 2 and pollutant emissions) are fully addressed or at least strongly mitigated. This purpose can be achieved by enhancing the efficiency of the combustion process and switching to the use of sustainable fuels. The most promising techniques thus far appear to be: Pre-Chamber (PC) combustion for Spark Ignition (SI) engines [9]; Dual Fuel (DF) combustion for Compression Ignition (CI) engines.
It is generally recognized by most researchers that DF combustion can provide ultralow NO x and soot emissions, and, at the same time, excellent Brake Thermal Efficiency (BTE) [10]. Moreover, low carbon-print fuels can be employed successfully in this combustion mode. In particular, Natural Gas (NG) [11], ethanol [12], methanol [13] and Liquefied Petroleum Gas (LPG) [14] can be used as low reactivity fuels, while biodiesel [15] and Dimethyl Ether (DME) [16] are suitable as high reactivity fuels.
DF NG-diesel combustion is particularly attractive since the hydrogen-to-carbon ratio of methane (CH 4 ) is the highest among hydrocarbons. The main drawback that affects this combustion technique is the high amount of UHC in the exhaust gas, especially at low loads. A viable strategy to address this issue is to add some H 2 to the lean premixed NG-air mixture entering the cylinder: even a small amount strongly promotes the oxidation of CO and UHC during combustion [17,18].
As an example, Rahnama et al. [19] explored the effect of H 2 addition to the intake air, varying the H 2 mole fraction from 0% to 5%, on a Heavy Duty (HD) DF NG-diesel engine. The analyzed operating condition is: 1300 rpm-Indicated Mean Effective Pressure (IMEP) equal to 9 bar (0.9 MPa). The results show a reduction of UHC emissions of about 67% (from 3.4 g/kWh to 1.1 g/kWh) and of CO emissions of about 46% (from 1.5 g/kWh to 0.8 g/kWh). Moreover, the authors report a reduction of soot of about 55%. The benefit of H 2 on soot emissions was confirmed by Zhou et al. [20], in similar conditions. Liu et al. [21] observed that H 2 bears a more noticeable effect on the early stage of the combustion process, when the higher combustion temperature associated with H 2 promotes the ignition of the surrounding NG.
Thripathi et al. [22] carried out a numerical study on the influence of the diesel injection strategy, considering a dual fuel, medium bore (102 mm) engine, running on hydrogen at 1500 rpm. In the single injection case, they found that advancing the injection timing reduces the HC, CO and soot emissions and increases the NO x emission and rate of pressure rise. Regarding dual injection, in order to match the targets of enhanced performance and reduced emissions, 16 • bTDC represented the optimum timing for pilot, 10% the optimum pilot mass, and 8 • bTDC the optimum main injection timing. The maximum value of gross Indicated Mean Effective Pressure (gIMEP) was about 10 bar (1.0 MPa) As can be deduced from the previous section, most of the published studies on DF combustion with hydrogen and/or NG are carried out on converted HD diesel engines running at 1300-1500 rpm. Conversely, the present paper is focused on a light-duty turbocharged compression ignition engine, running at a higher speed (3000 rpm) and comparable thermal load.
This type of engine may be of interest for a new generation of compact and efficient DF power units, whose application can vary from light commercial vehicles to industrial engines, delivering low-medium power rates. Standardization can reduce the cost of research and development, also permitting the transfer to advanced technologies already available for passenger cars (as an example, electronic injection systems, turbochargers, smart construction and manufacturing techniques . . . ). Moreover, in comparison to conventional SI engines running on NG, biogas, ethanol, hydrogen, et cetera, this new generation of light DF engines would bear the potential of a higher brake thermal efficiency, maintaining comparable dimensions and weight.
In order to convert such a potential into real value, it is necessary to gain a deep understanding of the DF combustion process, as well as to assess the influence of H 2 -NG blend composition. On the one hand, it would be highly desirable to run the engine on H 2 only (no CO 2 ); on the other hand, the issues of availability, storage and cost of H 2 in the short and medium term present a huge challenge for this proposition. Clearly, a trade-off must be identified, considering the evolving conditions. A practical proposition, at least for the foreseeable future, is to limit the volumetric share of H 2 to 30%: this figure represents the threshold to avoid major modifications to existing systems for gas distribution and injection.
In this context, the main goal of the current study is to provide an insight on DF combustion with NG and H 2 , considering light-duty diesel engines running at a representative speed and load. The influence of diesel injection strategy is assessed, along with the influence of the low reactivity fuel composition.
The work described in the paper consists of a 3D-CFD analysis, carried out by means of a customized version of KIVA-3V. As reported in a previous paper [23], the original diesel engine was modified in order to operate in DF NG-diesel mode. Then, a comprehensive experimental campaign was undertaken, at different operating conditions. The experimental data provided the basis for the validation of the 3D-CFD cylinder model. After validation, the model was employed to perform a set of simulations at an operating condition corresponding to 3000 rpm-Brake Mean Effective Pressure (BMEP) = 8 bar in DF H 2 /NG-diesel mode. The influence of an H 2 mole fraction spanning from 0% to 50% was investigated.

The Reference Engine
The engine used in the study is an automotive High-Speed Direct-Injection (HSDI) diesel unit, manufactured by FCA-VM Motori. It is a four cylinder in-line, with a total displacement of 2780 cc, four valves per cylinder. It is equipped with a high-pressure common-rail fuel-injection system, featuring six-holed injectors. The main characteristics of the engine are listed in Table 1. As mentioned above, the engine was modified in order to operate in DF NG-diesel mode: in detail, four NG injectors were installed on the engine inlet pipe, just downstream of the intercooler and about 500 mm before the intake manifold. The relatively long distance between the injectors and the intake manifold favors the formation of a homogeneous NG-air mixture and its uniform distribution among the 4 cylinders. Moreover, the EGR valve was constantly held shut, therefore the influence of charge dilution with exhaust gas was not considered. Different DF operating points, varying both the engine load and the amount of diesel substitution with NG, were investigated. The results of the experimental campaign are reported in [23].
The current study is focused on the operating condition corresponding to 3000 rpm BMEP = 8 bar. Table 2 reports the main parameters for both DF NG-diesel mode and for normal diesel mode (referred to as ND).

The 3D-CFD Model
For the combustion analysis, a customized version of the KIVA-3V code was used [24]. The code solves the conservation equations for evaporating fuel sprays, coupled with the 3-Dimensional Computational Fluid Dynamics (3D-CFD) equations of compressible, multicomponent, reactive gases in an engine cylinder with arbitrarily-shaped piston geometry. The DF combustion model implemented in the customized version of the code includes two different sub-models: the traditional Partially Stirred Reactor combustion mode (PaSR, fully described in [25]) and the Flame Propagation mode (fully described in [26]). The NG is modeled as a 3-component mixture of methane (CH 4 , 96 vol%), ethane (C 2 H 6 ) and propane (C 3 H 8 ), while the Diesel Oil Surrogate (DOS) model has been implemented in the code in order to represent the diesel fuel [27]. The DOS model is based on: a single component liquid fuel (equivalent chemical formula: H 14 H 28 ) containing the same main properties of real diesel and a two component vapor model constituted of a blend of n-heptane (H 7 H 16 ) and toluene (H 7 H 8 ). The mechanism developed to simulate the DF combustion is comprised of 81 species and 421 reactions.
The fuel models and the combustion mechanism were preliminary validated through comparison with experimental data in terms of ignition delay times in shock tube experiments and flame propagation data for constituent components of natural gas. More details and the validation process data are fully described in [26,28].
The customized version of the KIVA code includes also the Kelvin-Helmholtz/ Rayleigh-Taylor (KHRT) atomization model for the prediction of the diesel spray evolution during injection phase [29].
Finally, the same customized version of KIVA-3V code employed in the present study was already utilized by the authors to simulate different diesel and DF engines and the simulations results were consistently identified in good agreement with experiments, as reported in [27,28,[30][31][32].
The computational grid was built with the K3PREP pre-processor, included in the KIVA-3V software [24]. Due to the almost perfect axial-symmetric geometry of the combustion chamber and of the injection nozzle, characterized by six holes, a 60 • sector grid was considered and proper periodic boundary conditions were imposed.
In building the grid, particular attention was dedicated to assuring a good aspect ratio to the cells and to correctly reproduce the actual combustion chamber geometry. The squish region was adapted in order to precisely reproduce the engine compression ratio and a layer of six cells was imposed in the squish region at Top Dead Center (TDC).
In order to define the typical cell size, an analysis on different grids was carried out. Considering that the sector mesh is built rotating along the circumferential direction a 2D grid of the combustion chamber (see Figure 2c), the different grids were created varying both the reference cell dimension of the 2D mesh (from 0.5 to 1.7 mm) and the number of cells along the circumferential direction (from 15 to 35 cells, corresponding to an angular step between 1.7 • and 4 • per cell). According to K3PREP notation, in the following, the number of cells along the circumferential direction will be indicated as Ny. In total, 75 computation grids were built, with a number of cells at BDC varying from about 7.500 to about 188.000. Figure 2a,b summarizes the main results of the analysis: the former shows the dependency of one of the most important outputs of the simulations (the total heat released) as a function of typical cell size dimension for different Ny values; the latter illustrates the computational time for the different computational grids. Observing Figure 2a, it is possible to identify 2 different zones: for cell sizes greater than 1.1 mm, the output is dependent on the cell size itself, while, for cell sizes lower that 1.1 mm the results tend toward a single value (excluding only the points with Ny = 15, values that evidently represent an excessively coarse circumferential discretization). In particular, considering only the grids with Ny between 25 and 35, the results with a cell size lower than 0.8 mm are almost superimposed. Moreover, considering the computational cost, which increases exponentially with decreases in cell size (Figure 2b), it is possible to identify a cell size of 0.78 and Ny equal to 25 as the best tradeoff between accuracy and computational cost. This conclusion also aligns with previous analyses carried out by authors on other combustion geometries.
Ny equal to 25 as the best tradeoff between accuracy and computational cost. This conclu-sion also aligns with previous analyses carried out by authors on other combustion geometries.
The final computational grid, built with the criteria above explained, is reported in Figure 2c; it is comprised of about 80,000 cells at BDC and about 24,000 cells at TDC while the typical cell dimension of the 2D mesh varies between 0.5 and 1.0 mm (mean value 0.78 mm). The combustion simulations are carried out at closed valves, starting at Intake Valve Closing (IVC, −130 CA °AFTDC) and stopping at Exhaust Valve Opening (EVO, 112 CA °AFTDC). Initial conditions, such as pressure, temperature, trapped mass, and charge composition, are directly obtained from the experimental data. The initial flow field is imposed as a rigid vortex around the cylinder axis; its intensity (Swirl Ratio equal to 1.8) was calculated by a previous 3D-CFD intake stroke simulation. When assessing the influence of a diesel injection strategy, the restart function is used in order to minimize the computational effort: the compression stroke from IVC to −40 CA °AFTDC is simulated only once, while all the other calculations, characterized by different injection laws, start from the last point.  The final computational grid, built with the criteria above explained, is reported in Figure 2c; it is comprised of about 80,000 cells at BDC and about 24,000 cells at TDC while the typical cell dimension of the 2D mesh varies between 0.5 and 1.0 mm (mean value 0.78 mm).

Combustion Model Validation
The combustion simulations are carried out at closed valves, starting at Intake Valve Closing (IVC, −130 CA • AFTDC) and stopping at Exhaust Valve Opening (EVO, 112 CA • AFTDC). Initial conditions, such as pressure, temperature, trapped mass, and charge composition, are directly obtained from the experimental data. The initial flow field is imposed as a rigid vortex around the cylinder axis; its intensity (Swirl Ratio equal to 1.8) was calculated by a previous 3D-CFD intake stroke simulation. When assessing the influence of a diesel injection strategy, the restart function is used in order to minimize the computational effort: the compression stroke from IVC to −40 CA • AFTDC is simulated only once, while all the other calculations, characterized by different injection laws, start from the last point.

Combustion Model Validation
The first phase of the work consisted of the calibration and validation of the numerical model for both the ND and the DF operations at 3000 rpm-BMEP = 8 bar. The calibration basically consisted of tuning the shape of the diesel injection law, and some other minor parameters that were not defined in the set of initial conditions provided by the previous CFD-1D simulation. This parameter, referred to as gross Indicated Mean Effective Pressure (IMEP*), is calculated by the following formula: V is the engine unit displacement, p and V are pressure and volume within the cylinder.  Another parameter considered for the validation of the model is the specific work transferred from gas to piston in the interval from −40 CA • to 114 CA • after TDC (EVO). This parameter, referred to as gross Indicated Mean Effective Pressure (IMEP*), is calculated by the following formula: V d is the engine unit displacement, p and V are pressure and volume within the cylinder. Table 3 presents a comparison between simulation and experiments in terms of IMEP* and peak in-cylinder pressure. The data further confirm the consistency of the numerical model, in both combustion modes.

H 2 -NG Blend Simulations
In the second phase of the study, the calibrated numerical model was employed to assess the influence of composition of an H 2 -NG blend (also called Hythane) on performance and emissions at 3000 rpm/BMEP = 8 bar. As the H 2 volumetric fraction in the H 2 -NG blend varies, the global fuel energy contained in the premixed charge is kept constant. It is important to notice that, in most experimental studies on Hythane, the maximum volumetric concentration of H 2 is 30%. This limitation permits the replacement of Compressed Natural Gas (CNG) with Hythane with minimum modifications to the existing storage and injection systems. However, this numerical study was extended up to 50 vol% of H 2 , in order to provide a wider perspective.
The volumetric fractions of the various components used to represent the mixtures of Low Reactivity Fuels (LRF), are reviewed in Table 4. Further parameters are provided in Table 5, in particular: the mass of each main fuel component; the actual Air-Fuel (A/F) ratio of the premixed charge, (A/F) premix ; the correspondent stoichiometric value ((A/F) st, premix ) and the relative A/F ratio (λ premix ); the total A/F ratio (considering the injected diesel fuel), (A/F) tot ; the correspondent values of stoichiometric A/F ((A/F) st,tot ) and Lambda (λ tot ), where: Observing Table 5, one can notice that Lambda is much higher than 1, for both the premixed charge (1.83) and for the one including the injected diesel fuel (1.45). This means that the engine always works in lean conditions. Moreover, λ premix and λ tot are very similar for the different blends. This is due to the concurrent increment of (A/F) and (A/F) st as the fraction of H 2 increases in the H 2 -NG blend. The last two columns of Table 5 represent, respectively, the energy and the mass of H 2 divided by the energy and mass associated to the H 2 -NG-diesel mixture. Since the Lower Heating Value (LHV) of H 2 is more than twice that of NG, the mass fraction of H 2 is less than half its energy fraction.

Influence of Hydrogen Content and Start of Injection
In this section, the results of the combustion simulations with different H 2 volumetric fractions are presented and discussed. Figure 4a,b shows the in-cylinder pressure and the AHRR as a function of the volume fraction of H 2 in the H 2 -NG blend, ranging from 0% (baseline Dual Fuel (DFb) case) to 50%. As the amount of H 2 increases, the peak in-cylinder pressure always increases: the increment is about 12 bar (+20.5%) for the case with 30 vol% of H 2 (H30) and about 34 bar (+43%) for the case with 50 vol% of H 2 (H50). The increment of Peak Pressure Rise Rate (PPRR) illustrates another direct effect of the H 2 concentration growth: while PPRR for the DFb case is about 2.15 bar/CAD, this parameter rises at 3.40 bar/CAD for H30 and at 6.30 bar/CAD for H50 cas. Figure 4 also highlights the tight correlation between the peaks of AHRR and PPRR. As expected, the addition of H 2 accelerates the combustion process, shrinking its duration, without varying the Start Of Combustion (SOC) angle.

Influence of Hydrogen Content and Start of Injection
In this section, the results of the combustion simulations with different H2 volumetric fractions are presented and discussed. Figure 4a,b shows the in-cylinder pressure and the AHRR as a function of the volume fraction of H2 in the H2-NG blend, ranging from 0% (baseline Dual Fuel (DFb) case) to 50%. As the amount of H2 increases, the peak in-cylinder pressure always increases: the increment is about 12 bar (+20.5%) for the case with 30 vol% of H2 (H30) and about 34 bar (+43%) for the case with 50 vol% of H2 (H50). The increment of Peak Pressure Rise Rate (PPRR) illustrates another direct effect of the H2 concentration growth: while PPRR for the DFb case is about 2.15 bar/CAD, this parameter rises at 3.40 bar/CAD for H30 and at 6.30 bar/CAD for H50 cas. Figure 4 also highlights the tight correlation between the peaks of AHRR and PPRR. As expected, the addition of H2 accelerates the combustion process, shrinking its duration, without varying the Start Of Combustion (SOC) angle.  In order to limit the effects of H2 enrichment on peak in-cylinder pressure and PPRR, a sweep of the diesel Start Of Injection (SOI) angle was investigated. In detail, keeping the DFb case as a reference, the diesel injection law was advanced up to 3 CAD and delayed up to 6 CAD, by steps of 1 CAD. For each value of SOI, all H2-NG blends were simulated, for a total of 70 cases (10 SOI and 7 blends). The simulation results provided the basis for building the contour maps shown in Figures 5-8.
The colors of Figure 5a,b represent the IMEP* values as a function of H2 vol% and of diesel SOI. The solid red line indicates the value of IMEP* corresponding to the DFb case. In order to limit the effects of H 2 enrichment on peak in-cylinder pressure and PPRR, a sweep of the diesel Start Of Injection (SOI) angle was investigated. In detail, keeping the DFb case as a reference, the diesel injection law was advanced up to 3 CAD and delayed up to 6 CAD, by steps of 1 CAD. For each value of SOI, all H 2 -NG blends were simulated, for a total of 70 cases (10 SOI and 7 blends). The simulation results provided the basis for building the contour maps shown in Figures 5-8.
The colors of Figure 5a,b represent the IMEP* values as a function of H 2 vol% and of diesel SOI. The solid red line indicates the value of IMEP* corresponding to the DFb case.
The dashed lines of Figure 5a, superimposed on the IMEP* contour map, indicate the values of peak in-cylinder pressure, as a function of H 2 vol% and of diesel SOI; the yellow dashed line corresponds to the value of peak cylinder pressure of the DFb case. efficiency and at lower mechanical stress, in comparison to the DFb case.
Analyzing Figure 5b, with the PPRR curves (black dashed lines) plotted on the contour map, the curve corresponding to the DFb case is not shown as its value (2.15 bar/CAD) lies outside the domain of the graph. However, assuming a PPRR limit of 5 bar/CAD, most of the combinations H2 vol%-Delta SOI remain acceptable. It can be noticed that, accepting a peak cylinder pressure up to 110 bar, the limit of 5 bar/° can be reached with H30, at a much higher IMEP* (thus thermal efficiency) than DFb.
(a) (b) Figure 5. Contour maps of IMEP* for different H2 vol%-delta SOI combinations with peak in-cylinder pressure (a) and PPRR (b) levels superimposed.
Observing Figure 6a,b, where the trends of the combustion parameters CA50 and CA10-90 are plotted as a function of H2 vol%-Delta SOI, it is noted that combustion accelerates as the H2 mole fraction rises: at the SOI angle corresponding to DFb, both parameters decrease as H2 concentration increases. CA50 also strongly depends on diesel SOI: the earlier the start of combustion, primed by diesel injection, the smaller the CA50. Moreover, CA10-90 decreases as the SOI advances, since the higher combustion temperature favors the flame propagation (see Figure 6b).   Figure 7a,b displays BTE and combustion efficiency (η ), as a function of H2 vol% and Delta-SOI. The presence of H2 improves BTE. This is due to the combination of a faster combustion, associated with a higher efficiency of the thermodynamic cycle, and to a more complete combustion (higher η ). The last aspect is clearly visible in Figure 7b. Figure 8a-d shows the specific gaseous emissions, calculated as the ratio between the emissions mass flow rate and the engine brake power (P). P is evaluated as follows: where V is the total displacement in, n is the engine speed and η * is a correction parameter, calculated according to the following formula:   Figure 7a,b displays BTE and combustion efficiency (η ), as a function of H2 vol% and Delta-SOI. The presence of H2 improves BTE. This is due to the combination of a faster combustion, associated with a higher efficiency of the thermodynamic cycle, and to a more complete combustion (higher η ). The last aspect is clearly visible in Figure 7b. Figure 8a-d shows the specific gaseous emissions, calculated as the ratio between the emissions mass flow rate and the engine brake power (P). P is evaluated as follows: where V is the total displacement in, n is the engine speed and η * is a correction parameter, calculated according to the following formula: η * depends on the friction and pumping losses of the engine; it can be considered independent of the cases, for the reasons listed below: • friction losses mainly depend on engine speed, that does not change, and on in-cylinder peak pressures, that are quite similar; Observing Figure 5a, it is possible to notice that IMEP* rises as H 2 concentration and ignition advance increase. Keeping the same SOI of the DFb case, IMEP* raises from 11.7 bar to 13.0 bar (+11.3%) when passing from H0 to H30. A further increase of H 2 vol% to H50 causes IMEP* to reach up to 13.2 bar (+12.7%). On the other hand, the peak incylinder pressure and PPRR also increase, as previously discussed. A good trade-off is possible: moving along the red solid line (same IMEP* of DFb), H 2 can be increased up to 30% maintaining a slight benefit in terms of peak in-cylinder pressure; in a similar way, moving along the yellow dashed line (same peak cylinder pressure of DFb), H 2 can be increased up to 30% along with a slight increase of IMEP*.
Unfortunately, Figure 8d shows that the combustion process accelerates, as the H2 mole fraction increases: as a result, larger NOx emissions are discovered, due to the higher in-cylinder peak temperatures.
Moreover, SOI plays a fundamental role in pollutant emissions: an early start favors the combustion completion, reducing CO and UHC emissions, while increasing NOx, due to the higher temperatures.

Engine Performace and Emissions of Selected Operating Points
The contour maps of Figures 5-8 help to identify the most promising combinations of H2 mole fractions and diesel SOI. In particular, three operating conditions, all considering the maximum acceptable level of H2 in Hythane blends (30 vol%), were defined and compared to the DFb case.
The first operating condition was defined on the basis of Figure 5a: considering the "safety region" included between the red solid line and the yellow dashed line, a Delta_SOI of 3 CAD (delayed injection) is detected. This strategy permits the operation with a value of IMEP* higher than in the reference case, while complying with the limit of in-cylinder peak pressure of the reference case. As previously mentioned, since all combinations contain the same fuel energy as an input, an increase in terms of IMEP* describes an improvement of BTE. Such an operating point is referred to as "H30_SOI +3".
The second operating point is defined with consideration that PPRR is strongly linked to engine noise and vibrations [33]. In this study, the upper limit for PPRR is 5 bar/CAD, according to the literature. Therefore, the second operating point is detected at the intersection between the PPRR = 5 bar/CAD curve and the coordinate of 30 vol% H2. This corresponds to a further diesel injection advancement of 2 CAD and is highlighted by a red point in Figure 5b. The above-mentioned operating condition was referred to as "H30_SOI −2".
Finally, NOx emissions are considered: the third operating point, characterized by a 30 vol% H2 and Delta_SOI equal to 4.5 CAD (delayed injection), red point in Figure 8d Based on the previous considerations, the region between the red solid line and the yellow dashed line is identified as a "sweet spot", where the engine operates at higher efficiency and at lower mechanical stress, in comparison to the DFb case.
Analyzing Figure 5b, with the PPRR curves (black dashed lines) plotted on the contour map, the curve corresponding to the DFb case is not shown as its value (2.15 bar/CAD) lies outside the domain of the graph. However, assuming a PPRR limit of 5 bar/CAD, most of the combinations H 2 vol%-Delta SOI remain acceptable. It can be noticed that, accepting a peak cylinder pressure up to 110 bar, the limit of 5 bar/ • can be reached with H30, at a much higher IMEP* (thus thermal efficiency) than DFb.
Observing Figure 6a,b, where the trends of the combustion parameters CA50 and CA10-90 are plotted as a function of H 2 vol%-Delta SOI, it is noted that combustion accelerates as the H 2 mole fraction rises: at the SOI angle corresponding to DFb, both parameters decrease as H 2 concentration increases. CA50 also strongly depends on diesel SOI: the earlier the start of combustion, primed by diesel injection, the smaller the CA50. Moreover, CA10-90 decreases as the SOI advances, since the higher combustion temperature favors the flame propagation (see Figure 6b). Figure 7a,b displays BTE and combustion efficiency (η comb ), as a function of H 2 vol% and Delta-SOI. The presence of H 2 improves BTE. This is due to the combination of a faster combustion, associated with a higher efficiency of the thermodynamic cycle, and to a more complete combustion (higher η comb ). The last aspect is clearly visible in Figure 7b. Figure 8a-d shows the specific gaseous emissions, calculated as the ratio between the emissions mass flow rate and the engine brake power (P). P is evaluated as follows: where V d is the total displacement in, n is the engine speed and η * i is a correction parameter, calculated according to the following formula: η * i depends on the friction and pumping losses of the engine; it can be considered independent of the cases, for the reasons listed below: • friction losses mainly depend on engine speed, that does not change, and on in-cylinder peak pressures, that are quite similar; • pumping losses are related to the mass flow rate delivered by the engine, that is kept almost constant throughout the cases.
CO 2 specific emissions are expected to decrease as the H 2 mole fraction increases. The reasons include BTE improvement and the reduction of the carbon content in the H 2 -NG blend. This outcome is clearly visible in Figure 8a.
The completeness of combustion helps to also reduce CO and UHC (see Figure 8b,c). The result is coherent with Figure 7b, where η comb improves by increasing the H 2 content.
Unfortunately, Figure 8d shows that the combustion process accelerates, as the H 2 mole fraction increases: as a result, larger NO x emissions are discovered, due to the higher in-cylinder peak temperatures.
Moreover, SOI plays a fundamental role in pollutant emissions: an early start favors the combustion completion, reducing CO and UHC emissions, while increasing NO x , due to the higher temperatures.

Engine Performace and Emissions of Selected Operating Points
The contour maps of Figures 5-8 help to identify the most promising combinations of H 2 mole fractions and diesel SOI. In particular, three operating conditions, all considering the maximum acceptable level of H 2 in Hythane blends (30 vol%), were defined and compared to the DFb case.
The first operating condition was defined on the basis of Figure 5a: considering the "safety region" included between the red solid line and the yellow dashed line, a Delta_SOI of 3 CAD (delayed injection) is detected. This strategy permits the operation with a value of IMEP* higher than in the reference case, while complying with the limit of in-cylinder peak pressure of the reference case. As previously mentioned, since all combinations contain the same fuel energy as an input, an increase in terms of IMEP* describes an improvement of BTE. Such an operating point is referred to as "H30_SOI +3".
The second operating point is defined with consideration that PPRR is strongly linked to engine noise and vibrations [33]. In this study, the upper limit for PPRR is 5 bar/CAD, according to the literature. Therefore, the second operating point is detected at the intersection between the PPRR = 5 bar/CAD curve and the coordinate of 30 vol% H 2 . This corresponds to a further diesel injection advancement of 2 CAD and is highlighted by a red point in Figure 5b. The above-mentioned operating condition was referred to as "H30_SOI −2".
Finally, NO x emissions are considered: the third operating point, characterized by a 30 vol% H 2 and Delta_SOI equal to 4.5 CAD (delayed injection), red point in Figure 8d, leads to the same NO x specific emissions of the DFb case. The last point is referred to as H30_SOI +4.5.
The three operating points are compared in Figures 9-12. The DFb case is also reported and represented by a horizontal red line. Figure 9a-c displays IMEP*, η comb and BTE, respectively. As it can be seen, all the selected DF H 2 -NG/diesel cases exhibit better performances with respect to the DFb case. Moreover, as the diesel injection advances, both η comb and BTE improve and, as a consequence, the work output increases. Therefore, "H30_SOI −2" leads to a higher increment of IMEP* (11.7%), η comb (8.1%) and BTE (13.3%).
As expected, CA50 approaches closer to TDC as diesel injection advances (see Figure 10a). Furthermore, the combustion duration expressed by CA10-90 (Figure 10b), is always shorter than that of the DFb case, and decreases as the diesel injection advance increases. Therefore, the maximum variation from the DFb case in terms of CA10-90 corresponds to "H30_SOI −2" case (about −47%).
The three operating points are compared in Figures 9-12. The DFb case is also reported and represented by a horizontal red line. Figure 9a-c displays IMEP*, η and BTE, respectively. As it can be seen, all the selected DF H2-NG/diesel cases exhibit better performances with respect to the DFb case. Moreover, as the diesel injection advances, both η and BTE improve and, as a consequence, the work output increases. Therefore, "H30_SOI −2" leads to a higher increment of IMEP* (11.7%), η (8.1%) and BTE (13.3%).
As expected, CA50 approaches closer to TDC as diesel injection advances (see Figure  10a). Furthermore, the combustion duration expressed by CA10-90 (Figure 10b), is always shorter than that of the DFb case, and decreases as the diesel injection advance increases. Therefore, the maximum variation from the DFb case in terms of CA10-90 corresponds to "H30_SOI −2" case (about −47%). While the "H30_SOI +3" and "H30_SOI +4.5" cases contain the same, or even lower, peak in-cylinder pressure and PPRR compared to the DFb case, "H30_SOI −2" shows higher in-cylinder peak pressure and PPRR (see Figure 11a,b). In detail, comparing the of IMEP* (11.7%), η (8.1%) and BTE (13.3%).
As expected, CA50 approaches closer to TDC as diesel injection advances (see Figure  10a). Furthermore, the combustion duration expressed by CA10-90 (Figure 10b), is always shorter than that of the DFb case, and decreases as the diesel injection advance increases. Therefore, the maximum variation from the DFb case in terms of CA10-90 corresponds to "H30_SOI −2" case (about −47%). While the "H30_SOI +3" and "H30_SOI +4.5" cases contain the same, or even lower, peak in-cylinder pressure and PPRR compared to the DFb case, "H30_SOI −2" shows higher in-cylinder peak pressure and PPRR (see Figure 11a,b). In detail, comparing the Figure 10. Comparison between "H30_SOI −2", "H30_SOI +3", "H30_SOI +4.5" cases in terms of CA50 (a) and CA10-90 (b). Fuels 2022, 3, FOR PEER REVIEW 15 DFb case to the "H30_SOI −2" case, the in-cylinder peak pressure passes from 80 bar to about 110 bar, while PPRR increases by about 150% (from 2.2 bar/CAD to 5.1 bar/CAD).
Since the combustion process of the "H30_SOI −2" case is more advanced and faster, it is also characterized by higher NOx emissions with respect to the DFb case (+84%). Conversely, H30_SOI +3 and H30_SOI +4.5 show NOx emissions close to those of the DFb case.
Finally, all the selected DF H2-NG /diesel cases present lower CO, UHC and CO2 emissions, with respect to the DFb case. Due to the η and BTE trends, the maximum benefit is obtained advancing the diesel injection law of 2 CAD (H30_SOI −2). In this condition, CO and UHC tend toward zero (η is 100%), while CO2 emissions are reduced by about 12%, compared to the DFb case. Figure 11. Comparison between "H30_SOI −2", "H30_SOI +3", "H30_SOI +4.5" cases in terms of peak pressure (a) and PPRR (b). it is also characterized by higher NOx emissions with respect to the DFb case (+84%). Conversely, H30_SOI +3 and H30_SOI +4.5 show NOx emissions close to those of the DFb case.
Finally, all the selected DF H2-NG /diesel cases present lower CO, UHC and CO2 emissions, with respect to the DFb case. Due to the η and BTE trends, the maximum benefit is obtained advancing the diesel injection law of 2 CAD (H30_SOI −2). In this condition, CO and UHC tend toward zero (η is 100%), while CO2 emissions are reduced by about 12%, compared to the DFb case. Figure 12. Comparison between "H30_SOI −2", "H30_SOI +3", "H30_SOI +4.5" cases in terms of NOx (a), CO2 (b), CO (c) and HC (d).

Combustion Anlysis of Selected Operating Points
In order to gain a deeper insight into the influence of hydrogen addition and diesel SOI on DF combustion, the visualizations of O2 mass fraction and CH4 mass fraction in three cut planes, at various crank angles, are reported in the Appendix A and discussed in this section.
While the "H30_SOI +3" and "H30_SOI +4.5" cases contain the same, or even lower, peak in-cylinder pressure and PPRR compared to the DFb case, "H30_SOI −2" shows higher in-cylinder peak pressure and PPRR (see Figure 11a,b). In detail, comparing the DFb case to the "H30_SOI −2" case, the in-cylinder peak pressure passes from 80 bar to about 110 bar, while PPRR increases by about 150% (from 2.2 bar/CAD to 5.1 bar/CAD).
Since the combustion process of the "H30_SOI −2" case is more advanced and faster, it is also characterized by higher NO x emissions with respect to the DFb case (+84%). Conversely, H30_SOI +3 and H30_SOI +4. 5 show NO x emissions close to those of the DFb case.
Finally, all the selected DF H 2 -NG /diesel cases present lower CO, UHC and CO 2 emissions, with respect to the DFb case. Due to the η comb and BTE trends, the maximum benefit is obtained advancing the diesel injection law of 2 CAD (H30_SOI −2). In this condition, CO and UHC tend toward zero (η comb is 100%), while CO 2 emissions are reduced by about 12%, compared to the DFb case.

Combustion Anlysis of Selected Operating Points
In order to gain a deeper insight into the influence of hydrogen addition and diesel SOI on DF combustion, the visualizations of O 2 mass fraction and CH 4 mass fraction in three cut planes, at various crank angles, are reported in the Appendix A and discussed in this section.
In detail, Figure A1 reports the contour maps of O 2 mass fraction for "Diesel", "DFb" and "H30_SOI +3" cases in the following cut planes: • "O": cut plane normal to the cylinder axis; • "M": cut plane coincident with the symmetry plane of the sector mesh; • "B": cut plane coincident with the periodic boundary of the sector mesh.
As can be observed, the main difference among the diesel case and the DF cases ("DFb" and "H30_SOI +3 cases") is a higher and faster consumption of O 2 , due to the presence of the low reactivity fuel (NG or H 2 -NG blend, respectively) in the premixed charge. In particular, DF cases are able to exploit to a higher degree the O 2 present in the core and in the periphery of the combustion chamber compared to the diesel case.
Moreover, comparing "DFb" and "H30_SOI +3" cases, it appears the latter case is clearly characterized by a faster combustion process compared to the former, owing to the addition of hydrogen in the premixed charge. In fact, even if the SOI of the "H30_SOI +3" case is delayed by 3 CAD with respect to the "DFb" case, the consumption of O 2 in the "H30_SOI +3" case at 20 CAD ATDC is already higher compared to "DFb" case. Figure A2 illustrates the contour maps of CH 4 mass fractions for "DFb", "H30_SOI +3" and "H30_SOI −2" cases in the same cut planes depicted in Figure A1.
Observing Figure A2, the following considerations can be drawn. The oxidation of CH 4 in the diesel case occurs slower with respect to the DF cases. Furthermore, a higher fraction of CH 4 , located in the periphery of the combustion chamber, does not form part of the combustion process in the diesel case.
Finally, the higher the SOI advance, the faster (and more complete) the combustion process, as demonstrated by the comparison of "H30_SOI +3" and "H30_SOI −2" cases.

Conclusions
This paper explores the potential of H 2 addition in enhancing the combustion of a DF NG-diesel light duty engine. Distinct from other DF engines, based on Heavy Duty diesel units, this design is of particular interest, since it can be employed for a wider range of applications, from industry and agriculture to light commercial vehicles. In comparison to conventional SI engines running on NG, biogas, NG-H 2 blends, et cetera, it bears the potential for higher brake thermal efficiency.
A 3D-CFD model of the combustion system was constructed, calibrated and finally validated by comparison with experiments, then operated to investigate the effects of a partial substitution of NG with H 2 .
Simulations were carried out for different concentrations of H 2 in the H 2 -NG blend (up to 50 %vol), and by keeping the total amount of fuel energy provided by the premixed charge constant; the influence of the diesel injection strategy was also assessed (varying SOI).
Simulation results demonstrate H 2 enrichment accelerates the combustion process, improving its efficiency. As a consequence, BTE can increase by up to 12.6% for 50 %vol of H 2 and the same SOI of the baseline case, while specific UHC and CO emissions are dramatically reduced. Moreover, under the same operating conditions, CO 2 specific emissions are reduced by 19.7%, owing to the improvement of BTE and lower content of carbon in the premixed charge. The main drawbacks of hydrogen enrichment are: the increase of in-cylinder peak pressure, PPRR and NO x emissions.
However, a proper calibration of SOI permits the elimination of these shortcomings, maintaining an advantage, in comparison to standard DF NG-diesel combustion. As an example, considering a 30 %vol of H 2 and reducing the injection advance by 4.5 CAD, in-cylinder peak pressure, PPRR and NO x emissions remain unchanged, while BTE, UHC, CO, and CO 2 show improvement.    Figure A1. Comparison between "Diesel", "DFb" and "H30_SOI +3" cases in terms of O2 mass fraction visualized in different cut planes ("O": cut plane normal to the cylinder axis; "M": cut plane coincident with the symmetry plane of the sector mesh; "B": cut plane coincident with a cyclic boundary of the sector mesh), at different crank angles. Figure A1. Comparison between "Diesel", "DFb" and "H30_SOI +3" cases in terms of O 2 mass fraction visualized in different cut planes ("O": cut plane normal to the cylinder axis; "M": cut plane coincident with the symmetry plane of the sector mesh; "B": cut plane coincident with a cyclic boundary of the sector mesh), at different crank angles.   Figure A2. Comparison between "DFb", "H30_SOI +3" and "H30_SOI -2" cases in terms of CH4 mass fraction visualized in different cut planes ("O": cut plane normal to the cylinder axis; "M": cut plane coincident with the symmetry plane of the sector mesh; "B": cut plane coincident with a cyclic boundary of the sector mesh), at different crank angles.