Printer Friendly

Influence of Binary CNG Substitute Composition on the Prediction of Burn Rate, Engine Knock and Cycle-to-Cycle Variations.


Legal C[O.sub.2] emitting requirements and a worldwide increasing energy demand require a diversification on the fuel market, especially for automobile applications. When it comes to reaching the emission targets for passenger and freight transportation, natural and bio gases (CNG, Compressed Natural Gas and BNG, Bio Natural Gas) as well as synthetic methane based fuels (SNG, Synthetic Natural Gas) may play an important role. The advantages are well known compared to conventional fossil fuels: C[O.sub.2] savings of approximately 20 % can be realized compared to gasoline through a favorable H-to-C-ratio of methane. When adding renewable fuels, for example through electrolytically created hydrogen using excess solar power in daytime, the C[O.sub.2]-benefit can be even stronger and reach almost 100 % when using biogenically generated methane.

Further C[O.sub.2] saving potentials lie in the improvement of the CNG engine efficiency. On the one hand, the high knock resistance of methane allows to run an efficiency optimized 50 % mass fraction burnt point (MFB50) and higher compression ratios. On the other hand, dethrottling can be achieved by using higher EGR rates, which are possible due to the wide ignition limits of methane. To benefit from these advantages, the 0D/1D-simulation represents an important tool that serves reliable results with little effort, especially for transient operations.

Based on the necessary simplification of the reality with 1D-models, 1D-simulations heavily depend on the quality of their used submodels. For internal combustion engines it is of a high importance to model the combustion processes in great detail. Quasi-dimensional approaches, which are used to describe burn rates of natural-gas spark ignition engines, are mostly based on the modeling of the laminar flame speeds ([s.sub.L]). However, direct measurements of laminar flame speeds are usually performed for an air-fuel equivalence ratio range between [lambda] = 0.7 and 1.7 at pressures of only a few bar. Current approaches are extrapolated to unknown areas, which might cause contradictory data for laminar flame speed values.

To avoid problematic extrapolations for engine-related operation areas, reaction kinetics calculations are carried out to determine the laminar flame speed. The reaction mechanisms used for calculations follow known, physico-chemical principles which allow using a mechanism outside of its measurement-based validation range. Consequently, calculated laminar flame speeds can be approximated with a computing-time optimal correlation for 0D/1D-simulations.

When taking the influence of the fuel composition (for example rising amounts of hydrogen) on the laminar flame speed as well as knocking tendencies into account, the 0D/1D-simulation allows an even more detailed engine development process. In order to predict the operating range of lean running natural gas engines for example, determining knock and lean misfire limits is necessary. This becomes possible when additionally using a cycle-to-cycle variation model. Hence, emission reduction potentials and improvements in fuel consumption of different engine setups can be investigated.


The quasi-dimensional combustion model used in this study is presented in [1]. The laminar flame speed correlation proposed in this paper affects various equations regarding the combustion model. To clarify the influence, the principle idea of the combustion model will be outlined in the following.

Based on hemispheric flame propagation, the combustion chamber is divided into a burnt and an unburnt zone. Both zones are separated from each other by the flame front which is not considered as an additional, thermodynamical zone (see Figure 1).

The burn rate d[m.sub.b] is calculated by means of the mass [m.sub.F] entrained into the flame front and the characteristic burn-up time [[tau].sub.L] (Equation 1).

[[d.sub.mb]/dt] = [[m.sub.F]/[[tau].sub.L]] = [[m.sub.E]-[m.sub.b]/[[tau].sub.L]] (1)

The characteristic burn-up time depends on the laminar flame speed [s.sub.L] and the Taylor length [L.sub.T] (Equation 2).

[[tau].sub.L] = [[[tau].sub.L]/[s.sub.L]] (2)

The laminar flame speed for gasoline is calculated in accordance with [3]. For methane, the calculation is based on Gulder [4]. This calculation method extrapolates into unsurveyed ranges. Therefore, the laminar flame speed correlation proposed in this study substitutes the previous calculation method. Additionally, the influence of admixing ethane, propane, n-butane or hydrogen is implemented, enabling the calculation of binary CNG substitutes. The Taylor length can be computed with the global length scale l, the turbulence speed [u.sub.turb], the turbulent kinetic viscosity [v.sub.T] and the Taylor factor [[chi].sub.Taylor], with a value assumed to be 15 (Equation 3) [5].

[mathematical expression not reproducible] (3)

To calculate the turbulence speed [u.sub.turb], two different k-[member of] turbulence models can be used: a homogeneous, isotropic one as described in [5] and a quasi-dimensional model according to [6]. The latter model simplifies the tumble flow structure using a Taylor-Green eddy [7]. The resulting flow field yields the turbulence production using the same k-[member of] sub-models as are used in 3D-CFD but with simplified combustion chamber geometries. Thus, influences on the turbulent kinetic energy (TKE) such as changing cam durations and timings, as displayed in Figure 2, can be predicted. For combustion modeling, mainly the turbulence around top dead center firing (TDCF) is of great importance, where the model prediction matches the CFD simulation. Besides the influences on the Taylor length [l.sub.T] and the characteristic burn-up time [[tau].sub.L], respectively, [u.sub.turb] also affects the entrainment velocity [u.sub.E] (see Equation 4) which defines the entrainment mass flow d[m.sub.E] into the flame front (see Equation 5). As a consequence, the consideration of TKE-influencing measures, besides reliable laminar flame speed models, is of importance in order to properly predict the engine operating range.

[u.sub.E] = [u.sub.turb] + [s.sub.L] (4)

[d[m.sub.E]/dt] = [[rho].sub.ub] * [A.sub.fl] * [u.sub.E] (5)

The gas density in the unburnt zone ([[rho].sub.ub]) is based on the corresponding mass and volume. The flame surface [A.sub.fl] is calculated according to [1].


To evaluate the quality of the hitherto used Gilder-correlation [4], reaction kinetics calculations of laminar flame speeds are carried out for methane. Furthermore, [s.sub.L]-changes due to admixing ethane, propane, n-butane or hydrogen cannot be considered in the engine simulation yet. In this chapter, these secondary fuel components are going to be investigated so that influences can be understood and applied in quasi-dimensional burn rate models.

Reaction Kinetics Calculations and Results

Reaction mechanisms are the basis of reaction kinetics calculations. Such mechanisms are developed for specific fuels and contain the equations of all (known) elementary reactions taking place during the combustion. These reactions are studied in detail for a wide range of boundary conditions to investigate on their temperature and pressure dependency as well as on their material properties such as mass transfer coefficients of all molecules that are listed in the reaction equations.

To perform reaction kinetics calculations, a software is needed to use the information stored in a reaction mechanism. In this study, Cantera [9] is used and controlled via the computer language Python. Cantera offers one-dimensional flames as a calculation scenario, which is used to determine laminar flame speeds. Furthermore, the reaction mechanism GRI-Mech 3.0 [10] is included. This mechanism has been developed and validated to calculate natural gases and is widely used for methane calculations. However, it contains only proprietary information about propane and n-butane chemistry [10]. Hence, the performance of other mechanisms was tested to capture this admixture influence by comparing calculation results with measured sL values. For the admixture of propane, the USC C1-C3 mechanism [11] has been used. The USC C1-C4 [12] showed the best quality when calculating methane/n-butane mixtures.

It is important to note that with changing reaction mechanisms, the calculation results for pure methane diverge. As stated in [11], this divergence might result from the adaption of the USC mechanisms to represent higher hydrocarbons in addition to methane. To obtain consistent [s.sub.L]-values for pure methane, the USC mechanisms have only been used to calculate the relative influence of admixing propane or n-butane. This relative influence has been applied to the absolute [s.sub.L] values for pure methane, calculated with the GRI-Mech 3.0. Its results seem to be most trustworthy for pure methane in comparison with measurement data.

Via Python, it is possible to automate Cantera calculations over a wide range of boundary conditions. The boundary conditions used in this study are summarized in Table 1.

The different ranges have been chosen to cover all possible engine operation conditions. The values in the brackets were calculated, but not again used, when applying the model approach to the results of the reaction kinetics calculations. They served as a validation basis for the model behavior outside of its adaption range. The blending rate limit of 40 mol-% results on the one hand from the non-linear behavior of [H.sub.2]-admixture above 40 mol-%, where the transition from methane-controlled to hydrogen-controlled combustion begins [13]. On the other hand, a mixture of 60 mol-% methane and 40 mol-% hydrogen represents a methane number of 60, which already indicates a very high knocking tendency. Besides this, mixtures with considerably low methane numbers could be calculated by using high amounts of n-butane, for example 40 mol-% is equal to a methane number of 32.1.

Figure 3 shows the calculated laminar flame speeds of methane for varying residual exhaust gas fractions at different temperatures and pressures. While higher temperatures increase the speed of chemical reactions and, thus, the laminar flame speed, higher pressures or residual exhaust gas fractions show an inhibiting influence. Besides increasing the heat capacity, and thereby lowering the temperature, the residual exhaust gas dilutes the combustible mixture. This dilution decreases the fuel and the oxidizer concentration and, therefore, impedes the chemical reactions.

The pressure influence on the flame speed can be explained by Le Chatelier's principle [14]: High pressures promote reactions with fewer products than reactants. With a lower number of molecules, the pressure is reduced. For combustion reactions, this can be translated to the acceleration of chain-breaking reactions. In these reactions, fuel radicals recombine to a single stable molecule, hence breaking the chain reaction of the combustion and consequently decreasing the flame speed.

The calculated influence of rising hydrogen amounts for different fuel-air equivalence ratios [PHI] (=1/[lambda]) is compared with measurement data from [15] in Figure 4. With a higher amount of [H.sub.2], the laminar flame speed of the mixture increases significantly due to the high reactivity of hydrogen molecules and matches the measured trend. For mixtures of methane with ethane, propane or n-butane, the same qualitative trends of calculation results can be observed. However, the absolute influence on the flame speed is less significant, which matches with the trends measured in [16]. On the contrary, the knocking tendency is strongly increased by adding even small amounts of n-butane to methane for example, as described in the later chapter "Fuel Composition Influence on Knock Simulation".

When varying [PHI], a maximum in flame speed for slightly rich mixtures can be observed. This is a result of a chemical equilibrium, which only allows a conversion of all available oxygen with a surplus of fuel. With leaner or richer mixtures, the concentration of either reactant or oxidant decreases, which results in a lower speed of the combustion reactions. Additionally, the excessive amount of air or fuel has to be heated up, which brings a lower temperature and, thus, reduces the laminar flame speed.

It has to be noted that among other effects, laminar flame speed measurements might get influenced by flame wrinkling, for which a mathematical correction is necessary, eventually causing measurement uncertainties. Furthermore, flame wrinkling effects limit the range of available measurements to relatively low temperatures and pressures.

Figure 5 compares the influences of residual gas (EGR) and [lambda] on [s.sub.L] for pure methane at a pressure of 1 bar. The x-axes are scaled to match the curves of residual gas and [lambda] for 300 K. The difference in the scale indicates the different influences of [lambda] and EGR, since [lambda] = 1.3 and a residual gas fraction of 22.1 mass-% represent the same degree of charge dilution. On the one hand, this results from different heat capacities of the air and the residual gas, which influence the heating of the unburnt gas as well as the flame temperature. On the other hand, the different chemical behavior of the reactive excessive air at [lambda] = 1.3 changes the reactions that take place during the combustion, compared to the nearly inert residual gas. The increasing difference between the blue and the red lines with higher temperatures represents the temperature dependency between [lambda] and the residual gas influence, which needs to be taken into consideration when approximating the correlation for the laminar flame speeds.

In general, the calculated laminar flame speeds of all investigated mixtures show plausible trends over varied boundary conditions, matching available measurement data inside the range of measurement uncertainty. Since the single reactions that build a reaction mechanism can be studied in a far bigger boundary condition range than the laminar flame speeds, it is possible to use reaction mechanisms outside of their validation range [17]. Due to the lack of measurements, this range is relatively small. Therefore, only results of reaction kinetics calculation allow creating a laminar flame speed correlation for the engine application.

Model Approach

Figure 6 compares the laminar flame speeds that result from the Gilder-correlation [4], as described by the Equations 6 and 7, with reaction kinetics calculation results for pure methane at [lambda] = 1.

[mathematical expression not reproducible] (6)

[S.sub.L,0] ([PHI]) = 0.422 * [[PHI].sup.0,15] * [e.sup.-5.18]*[([PHI]-1.075).sup.2] (7)

Although the basic trends of pressure and temperature are similar, the absolute values from the Gilder-correlation are significantly lower, particularly at high temperatures. These boundary conditions are very important for burn rate calculations, as the combustion sets in at relatively high temperatures. A decreasing pressure results in a stronger increase of [s.sub.L] compared to the calculation results.

Besides the direct influence of temperature, pressure, air-fuel equivalence ratio, residual exhaust gas fraction and fuel composition on the laminar flame speed, further cross influences of different boundary conditions can be observed. The most prominent influence is that of the temperature on the pressure or vice versa. Investigations according to [18] propose a pressure dependency of the parameter [alpha] and a temperature dependency of [beta] (Equation 6) According to [18], these parameters also should be [lambda]-dependent, whereby the [lambda]-influence shows a non-derivable trend. Additionally, Equations 6 and 7 make it impossible to reproduce temperature trends that are shown in Figure 6. This indicates that the Gilder-correlation requires a revision in order to account for these influences. As a consequence, the parameters of another correlation for laminar flame speeds of iso-octane, described in [19] are adapted to match the calculation results of [s.sub.L] for methane. The correlation considers reaction kinetics as well as several cross influences, and is very promising in reaching a higher accordance between the reaction kinetics calculations and the correlated values. The correlation is described by the Equations 8 to 13.

[mathematical expression not reproducible] (8)

A([T.sup.0]) = F * [e.sup.-G / [T.sup.0]] (9)

[mathematical expression not reproducible] (10)

[mathematical expression not reproducible] (11)

[mathematical expression not reproducible] (12)

[mathematical expression not reproducible] (13)

All parameters that are not described by any sub-equations are calibration parameters. To match reaction kinetics calculation results for pure methane, the parameters have been calibrated by the least error square method for the boundary conditions that are listed in Table 1. To achieve the best possible conformity, errors at boundary conditions that are relevant for the engine operation have been weighted stronger than the not relevant ones, for example those at high pressures and at 300 K. The resulting parameters are listed in Table 2 and 3. [S.sub.1] to [S.sub.4] are splines, depending on Z*. During the calibration process, Equation 13 was changed by adding the exponent c to better match the influence of EGR on the calculated burnt temperature [T.sub.b].

The adaption quality (deviation in percent) between the correlation and the calculation for pure methane is displayed in Figure 7 at [lambda] = 1.

There are two noticeable areas of high deviations. At very high pressures above 170 bar and relatively low temperatures, a difference of 18 % and higher appears. For a combustion engine, however, the combination of these pressures and temperatures never takes place. For higher temperatures, the difference decreases. It has to be noted that the correlation adaption has only been carried out for pressures up to 100 bar. Higher pressures, hence, represent an extrapolation area, that shows a fairly well match between the calculation and the correlation. The same applies for the air-fuel equivalence ratio [lambda], where a good conformity between the correlation and the calculation can be observed throughout to the flammability limits although the correlation only had been adapted for [lambda] [less than or equal to] 1.7.

The big difference at pressures around 25 bar and low temperatures is a result of compromises that had to be made during the adaption of the correlation parameters. It is possible to reach a better conformity here, but other areas would, by that, be affected adversely. The white solid line represents the temperature and pressure trace of the unburnt zone (compare chapter "Burn Rate Modeling") during the high pressure cycle at full load, beginning with "cycle start" and ending with "cycle end". The white arrows mark the start of combustion as well as the MFB95 point, where 95 % of the fuel has already burnt. It can be seen that the main part of the combustion lies inside areas of low deviation between the calculation and the simulation, which explains the validity of the compromises that were made. When reducing the load, the white curve shifts towards lower pressures with similar temperatures. The displayed temperature/pressure trace, hence, represents the worst case. When varying [lambda] or the residual exhaust gas content, the level of deviation remains similar; most T/p-combinations stay below 12 %.

The trends displayed in Figures 6 and 7 can also be found in Figure 8, showing the calculated laminar flame speed, the Gilder-correlation and the correlation proposed in this work for one engine cycle. With Gilder strongly differing from the reaction kinetics calculation results, the general trend for [s.sub.L] can be reproduced much better when using the new correlation. Compared to the calculation results, the maximum error is approximately 14 %.

Since the correlation proposed in [19] only accounts for pure fuels, an extension of the formulas is necessary to include the influence of ethane, propane, n-butane and hydrogen admixtures. For the alkanes, a dependency of the calibration parameters F, [n.sub.a] and [n.sub.EGR] on the secondary component mass fraction [??] as described in Equations 14 to 16, is sufficient to find a good match between the correlation and the calculated results.

[mathematical expression not reproducible] (14)

[mathematical expression not reproducible] (15)

[mathematical expression not reproducible] (16)

The used calibration parameters are listed in Table 4 dependent on the alkane type. Just as with pure methane, they were determined by minimizing their weighted error squares.

The admixture-dependent parameters with the index [C.sub.x][H.sub.y] defined in the Equations 14 to 16, replace the existing parameters F, [n.sub.a] and [n.sub.EGR] in Equations 9 to 11. If [??] is zero, the equations are reduced to the base calibration parameters for pure methane, as listed in Table 2. Since the splines [S.sub.1] to [S.sub.4] are independent from the content of the secondary component, the equation for Z*, given in [19], is redefined with Equation 9, where [L.sub.min,CH4] is used. Hence, it can be avoided to change any Z*-values when changing the fuel composition.

The resulting difference in percent between the reaction kinetics calculation and the expanded correlation for 40 mol-% alkane admixture is shown in the Figures 9 to 11. In general, the trends and the degrees of the deviation are comparable to pure methane, which underlines the high quality of the expanded correlation.

For hydrogen, it is necessary to take the pressure- and [lambda]-dependency of hydrogen admixture into consideration. For this idea, Equation 17 was developed, which uses the calibration parameters that are listed in Table 5. Similar to the correlation that was developed for alkane admixtures, Equation 17 is reduced to the parameter F for pure methane if [??] is zero.

[mathematical expression not reproducible] (17)

Figure 12 shows the deviation between the calculated results and the extended correlation for a mixture between 60 mol-% of C[H.sub.4] and 40 mol-% of [H.sub.2]. When comparing the results with Figure 7, only a slight change in deviation distribution and no change in deviation level become apparent which underlines the good quality for the extension of the formula. In contrast to this, a simple adaption of existing correlation parameters with changing fuel composition, as it is sufficient for alkane admixture, would result in an increase of the deviation level.

Test in the Burn Rate Model

To validate the updated and expanded correlation not only in terms of conformity with calculated flame speeds, but also with its influence on burn rate simulations, it is implemented in the burn rate model, which was described in the chapter "Burn Rate Modeling". For this purpose, measured pressure traces of a single-cylinder research engine with 0.6 l of displacement volume for methane and methane with 30 mol-% hydrogen at 2000 rpm and full load are analyzed in order to calculate the burn rate by performing a pressure trace analysis (PTA). These burn rates are compared to simulation results in Figure 13, where the model is calibrated for pure methane.

While the shapes of simulated burn rates slightly differ from the measurements, the general influence of hydrogen admixture is well predicted: With a higher laminar flame speed, the burn rate increases faster and reaches higher peak values. As a result, the burn duration decreases, since the change in mass-specific heating value when admixing hydrogen is considered to reach fuels with a comparable energy content and, thus, resulting in similar integrated burn rates.

When using the Gilder correlation with the same model calibration parameters and 50 % mass fraction burnt point, the burn rate increases slower and reaches lower peak pressures due to lower laminar flame speeds. Consequently, the burn duration increases. At present, this problem could get solved by increasing the turbulence level.


The simulation of engine knock is based on an empirical approach, considering the complex process of auto-ignition in the unburnt zone as a collective reaction. This is in favor for the computation time compared to kinetic model approaches. The temperature in the unburnt zone is assumed to be homogenous and the global reaction rate as well as the ignition delay [tau] can be described by the Arrhenius equation [20]:

[tau] [approximately equal to] c * [p.sup. -a] * [[E.sub.A]/eRT] (18)

For the knock simulation, the parameters a, c and the activation energy [E.sub.A] must be determined as a reference point at a predefined knock rate limit and the Arrhenius equation has to be transferred to the combustion process of the engine and its temperature and pressure history. Therefore, the formulation of Franzke [21] for the knock integral value [I.sub.k] is used, which can be considered as the energy level of the unburnt gas. As [(x).sub.c] is the critical concentration of radicals causing auto-ignition, the [I.sub.k] describes the integral at which the state of knock is reached [21] according to the Equation:

[mathematical expression not reproducible] (19)

The investigation of fuel composition influences on knock simulation is modeled according to the FVV project [22]. The analysis of the single cycle pressure signal described in [22] shows a very early knock onset of MFB40 to MFB60 which stands in contrast to the current knowledge about Otto fuel knock onsets as they are usually expected at MFB75 to MFB90 [20]. Due to the chemical reaction inertia of the compact C[H.sub.4] molecule [23], it is supposed that there is no prompt conversion of the overall unburnt mass. A spontaneous reaction of more than 50 % of the unburnt mass would lead to very high pressure gradients and cause massive engine damage.

In order to calibrate the knock model in a first step, a knocking operation point of 100 % methane is evaluated, which is here considered as a reference gas. The basic parameters are subsequently applied to the binary gas mixtures. As a result, the knock integral values at 5 % knock ratio limit differ significantly, depending on the secondary gas and the methane number. This can be explained by the fuel impact, which is represented in the Arrhenius equation with the activation energy [E.sub.A]. To receive a uniform [I.sub.K] value at the knock limit, the [E.sub.A] values are calculated in dependency of the methane number as shown in Figure 14. The remaining parameters of Equation 18 are set to a = -1.1 and c = 5.58, R represent the universal gas constant. For a fixed methane number, the C[H.sub.4]/[H.sub.2] blends reach the highest activation energies that translate to the lowest knock tendencies. Considering the methane-alkane blends, activation energies decline with increasing chain length, which corresponds with the experimental results. The correlation between methane number and the activation energy can be described as a linear polynomial for all binary gas mixtures.


Compared to compression ignition engines, the modeling of sparkignition engine burn rates is much easier due to relatively homogeneous mixture distribution. Using entrainment-approaches (compare chapter "Burn Rate Modeling"), high prediction qualities can be achieved, which results in a good reproducibility when influences such as residual gas, charge dilution, charge motion and turbulence due to variable valve trains are applied.

Models that are capable of accounting for these influences allow simulating for example fuel consumption-optimal valve timings for the lean operating limit of stationary engines. To check if the application of the determined valve timings to a real engine is possible, it is important to take cycle-to-cycle variations (CCV) into consideration. This can be done by using the CCV model that is developed in [24]. It is based on the burn rate model, described in chapter "Burn Rate Modeling" and it applies stochastic noises of a constant bandwidth to the turbulence and the mixture distribution. Depending on boundary conditions such as charge dilution, engine speed, turbulence and temperature, the applied noises cause fluctuations in sub-models, which result in fluctuations of the indicated mean effective pressure (IMEP). This allows predicting statistical IMEP fluctuations, represented by the coefficient of variation (COV IMEP). If the COV IMEP exceeds a certain limit, the CCV are considered too high for a stable engine operation.

Figure 15 shows a comparison of measured and simulated cycle-to-cycle variations of an MTU BR4000 natural gas engine (57.2 l of displacement volume) for air fuel ratio or combined air fuel ratio/ignition point variations. It is apparent here that the effects of leaning are well reproduced by the CCV model. In addition, the stabilizing effect of the ignition timing on the fluctuation level is modeled satisfactorily. This proves the applicability of the CCV model to be independent on engine size or on fuel type.


The combination of all models described in the previous chapters allows predicting the operating ranges for spark-ignition engines, accounting for the influence of knock on the one hand and unstable engine operation due to cycle-to-cycle variations on the other hand.

For stationary CNG engines, the prediction of the operating range is of special interest in terms of in-cylinder nitric oxide (NO) emission reduction by lean engine operation and high efficiency through earliest possible MFB50. The knock limit restricts the operating range for an early MFB50 and a lean misfire limit restricts for a late MFB50. These limits for pure methane are represented in Figure 16 by the blue lines. When exemplarily adding 30 mol-% hydrogen (which equals a methane number of 70), the burn duration would decrease. If the fuel influence on knock is neglected (compare chapter "Fuel Composition Influence on Knock Simulation"), the decreased burn duration allows setting an earlier MFB50, which is represented by the dashed red line on the left. This qualitative behavior also applies for the other investigated fuel mixtures that have a methane number of 70. Since the quantitative difference between the fuels is small, their individual knock limits with unchanged activation energy [E.sub.A] are not displayed.

When additionally regarding composition influences on knock by changing [E.sub.A], the knock limit is shifted towards late MFB50, which is represented by the red dashed line in the middle for a C[H.sub.4]/[H.sub.2] blend. This shows that the influence of hydrogen admixture on knock overcomes the positive influence on burn duration, which in turn narrows the engine operating range. Concerning alkane admixtures, this effect is even more prominent, which is consistent with the trend displayed in Figure 14. In contrary, the stabilizing effect of higher laminar flame speeds on the combustion allows higher mixture dilutions, which widens the operating range. This behavior is represented by the shifted lean misfire limits towards late MFB50, illustrated in Figure 16. As underlined by the close-up from CAD 205 to 215, the degree of shift correlates with the level of laminar flame speed. Hydrogen or ethane admixtures cause a similarly high rise in laminar flame speed and, thus, allowing later MFB50 or higher [lambda]-values than [C.sub.3][H.sub.8] / C[H.sub.4] or n-[C.sub.3][H.sub.10] / C[H.sub.4] blends, respectively, whose flame speeds are similar but lower compared to C[H.sub.4] / [H.sub.2]. Minor differences in predicted lean misfire or knock limits might on the one hand result from [s.sub.L]-correlation imperfections and on the other hand from changing boundary conditions. That is caused by changing in-cylinder masses at the inlet valve closing event due to different air requirements for different fuels when keeping the fuel energy constant at constant [lambda].

Although the influence of secondary component admixture on lean misfire limit seems relatively small concerning the MFB50, the possibility to significantly increase [lambda] can result in decreasing NOx-emissions. For an individual application, this benefit has to counterbalance the reduced operating range due to shifted knock limits. In this context, hydrogen represents the best compromise.


Based on reaction kinetics calculation, a new correlation for laminar flame speeds of methane is proposed in this paper. The correlation aims to replace the Gilder correlation, which is still in use today for engine simulations. The proposed correlation is expanded to account for secondary fuel components such as ethane, propane, n-butane or hydrogen. This expansion allows investigating the influence of binary methane-based CNG substitutes on the combustion. In combination with the presented knock model adaption, changes in the fuel compositions, which for example differ among multiple natural gas sources, can be implemented for a simulation. Hence, the minimum fuel quality in terms of knock can be determined. When additionally using a cycle-to-cycle variation model, the full engine operating range for varying fuel compositions can be calculated. On the one hand, this allows evaluating available engine control ranges for different [lambda]-values. On the other hand, raw NO emission reductions, for instance, can be located.

In general, the presented study shows how the collaboration of several 0D/1D models and their mutual influences allow to predict engine operation limiting factors, which improves a computer aided engine development process significantly. Particularly for stationary gas engines, where experimental investigations are often expensive and limited to single-cylinder research engines, 0D/1D simulations offer a high potential in reducing monetary effort and speed up the development process by considering the behavior of a full engine.


[1.] Grill, M., Billinger, T., and Bargende, M., "Quasi-Dimensional Modeling of Spark Ignition Engine Combustion with Variable Valve Train," SAE Technical Paper 2006-01-1107, 2006, doi:10.4271/2006-01-1107.

[2.] Wenig, M., Grill, M., and Bargende, M., "A New Approach for Modeling Cycle-to-Cycle Variations within the Framework of a Real Working-Process Simulation," SAE Int. J. Engines 6(2):1099-1115, 2013, doi:10.4271/2013-01-1315.

[3.] Heywood, J. B., "Internal Combustion Engine Fundamentals,"McGraw-Hill Series in Mechanical Engineering, 1988.

[4.] Gulder, O., "Correlations of Laminar Combustion Data for Alternative S.I. Engine Fuels," SAE Technical Paper 841000, 1984, doi:10.4271/841000.

[5.] Bargende, M., "Ein Gleichungsansatz zur Berechnung der instationaren Wandwarmeverluste im Hochdruckteil von Ottomotoren," Ph.D. thesis, TH Darmstadt, 1991.

[6.] Bossung, C., "Turbulenzmodellierung fur quasidimensionale Arbeitsprozessrechnung," final report for FVV-project no. 1067, Frankfurt am Main: Forschungsvereinigung Verbrennungskraftmaschinen, 2014

[7.] Taylor, G. I. and Green, A. E., "Mechanisms of the Production of Small Eddies from Large Ones," Proceedings of the Royal Society of London. Series A, Mathematical and Physical Sciences, Vol. 158, No. 895, 1937

[8.] Bossung, C., "Quasidimensionales Ladungsbewegungs- und Turbulenzmodell fur die Motorprozessrechnung," Presented at ATZ live - Ladungswechsel im Verbrennungsmotor, 2014

[9.] Goodwin, D. G., Moffat, H. K., and Speth, R.L., "Cantera: An Object-oriented Software Toolkit for Chemical Kinetics, Thermodynamics, and Transport Processes,", 2015. - Version 2.2.0

[10.] Smith, G. P., Golden, D. M.; Frenklach, M., Eiteener, B.; GRI-Mech 3.0. 2000

[11.] Qin, Z., Lissianski, V. V., Yang, H., Gardiner, J., "An Optimized Reaction Model of C1-C3 Combustion," Version: 2000

[12.] Wang, H., You, X., Joshi, A. V., Davis, S. G. , "USC Mech Version II. High-Temperature Combustion Reaction Model of H2/CO/C1-C4 Compounds," Version: May 2007

[13.] Di Sarli, V. and Di Benedetto, A., "Laminar burning velocity of hydrogen-methane/air premixed flames," International Journal of Hydrogen Energy 32(2007), Nr. 5, S. 637-646

[14.] Turns, S. R., "An introduction to combustion: concepts and applications, 2nd Edition"(McGraw-Hill series in mechanical engineering, 2000), ISBN 0-07-230096-5.

[15.] Hermanns, R., "Laminar Burning Velocities of Methan-Hydrogen-Air Mixtures," Ph.D. thesis, Technische Universiteit Eindhoven, 2007

[16.] Dirrenberger, P., Le Gall, H., Bounaceur, R., Herbinet, O. et al, "Measurements of Laminar Flame Velocity for Components of Natural Gas," Energy and Fuels 25 (9): 3875-3884, 2011

[17.] Warnatz, J., Maas, U. and Dibble, R. W., "Verbrennung - Physikalisch-Chemische Grundlagen, Modellierung und Simulation, Experimente, Schadstoffentstehung", 3. Edition,Springer-Verlag Berlin Heidelberg, 2001

[18.] Konnov, A., "The Temperature and Pressure Dependences of the Laminar Burning Velocity: Experiments and Modelling" In: Proceedings of the European Combustion Meeting - 2015. Budapest, Hungary, 2015

[19.] Ewald, J., "A Level Set Based Flamelet Model for the Prediction of Combustion in Homogeneous Charge and Direct Injection Spark Ignition Engines," Ph.D. thesis, Rheinisch- westfalische Technische Hochschule Aachen, 2006

[20.] Spicher, U. and Worret, R., "Entwicklung eines Kriteriums zur Vorausberechnung der Klopfgrenze," FVV-Abschlu[ss]bericht, Vorhaben Nr. 700, 2002.

[21.] Franzke, D., "Beitrag zur Ermittlung eines Klopfkriteriums der ottomotorischen Verbrennung und zur Vorausberechnung der Klopfgrenze," Ph.D. thesis, TU Munchen, 1991.

[22.] Scharlipp, S. and Urban, L, "Methan-Kraftstoffe: Potenzialstudie und Kennzahlen," final report for FVV-project no. 1126, Frankfurt am Main: Forschungsvereinigung Verbrennungskraftmaschinen, 2015

[23.] Van Basshuysen, R., "Erdgas und erneuerbares Methan fur den Fahrzeugantrieb. Wege zur klimaneutralen Mobilitat," 2015.

[24.] Wenig, M., "Simulation der ottomotorischen Zyklenschwankungen, "final report for FVV-project no. 995, Frankfurt am Main: Forschungsvereinigung Verbrennungskraftmaschinen, 2012


FKFS - Research Institute of Automotive Engineering and Vehicle

Engines Stuttgart

Pfaffenwaldring 12

D-70569 Stuttgart, Germany

IVK - Institute of Internal Combustion Engines and Automotive


Pfaffenwaldring 12

D-70569 Stuttgart, Germany


M.Eng. Sebastian HANN

Dr.-Ing. Michael GRILL

Dipl.-Ing. Lukas URBAN

Prof. Dr.-Ing. Michael BARGENDE


The sub-models from [24], [6] and [22] are the result of research tasks defined by the Forschungsvereinigung

Verbrennungskraftmaschinen e.V. (FVV, Frankfurt) and conducted at the Institut fur Verbrennungsmotoren und Kraftfahrwesen (IVK) of the University of Stuttgart. They were self-financed by FVV ([24]) or within the framework of a program to promote cooperative industrial research (Industrielle Gemeinschaftsforschung und -entwicklung, IGF) by the Bundesministerium filr Wirtschaft und Technologie (BMWi) of the Federal Republic of Germany by means of the Arbeitsgemeinschaft industrieller Forschungsvereinigungen e. V. (AiF), IGF-Nr. 17143 N ([6]) and IGF-Nr. 17573 N ([22]), based on a decision of the German Federal Parliament.

Working groups accompanied the research works of [24], [6] and [22]. The authors would like to thank these working groups and the companies involved for their support and the BMWi, AiF and FVV for providing financing.


[A.sub.fl] - Flame surface

BMEP - Break mean eff. pressure

[[chi].sub.Taylor] - Taylor factor

C - Carbon atom

[C.sub.2][H.sub.6] - Ethane

[C.sub.3][H.sub.8] - Propane

CCV - Cycle-to-cycle variation

CFD - Computational fluid dynamics

C[H.sub.4] - Methane

CNG - Compressed natural gas

C[O.sub.2] - Carbon dioxide.

COV - Coefficient of variation

[dm.sub.b] - Mass entering burnt zone

[dm.sub.E] - Mass entering flame front

[epsilon] - Turbulence dissipation

[E.sub.A] - Activation Energy

EGR - Residual exhaust gas

H - Hydrogen atom

[H.sub.2] - Hydrogen

[I.sub.k] - Knock integral

IMEP - Indicated mean eff. pressure

[lambda] - Air-Fuel equivalence ratio

l - Global length scale

[l.sub.T] - Taylor length

[m.sub.F] - Mass in flame front

MFB50 - 50 % mass fraction burnt

n-[C.sub.4][H.sub.10] - n-Butane

NO - Nitric oxide

[[nu].sub.T] - Turbulent kinetic viscosity

[PHI] - Fuel-air equivalence ratio

[p.sub.0] - Reference pressure

PTA - Pressure trace analysis

[[rho].sub.ub] - Density of unburnt zone

[s.sub.L] - Laminar flame speed

[tau] - Ignition delay

[T.sub.0] - Reference temperature

[T.sup.0] - Inner layer temperature

[T.sub.b] - Temperature burnt zone

[[tau].sub.L] - Characteristic burn-up time

TDCF - Top dead center firing

TKE - Turbulent kinetic energy

[T.sub.u] - Temperature unburnt zone

[u.sub.E] - Entrainment velocity

[u.sub.turb] - Turbulence speed

X - Mole fraction

[(x).sub.c] - Critical radical concentration

Y - Mass fraction

[Z*] - Z* for [lambda] = 1

Sebastian Hann, Lukas Urban, and Michael Grill


Michael Bargende

IVK, University of Stuttgart

Table 1. Boundary conditions for reaction kinetics calculations

(Unburnt) Temperature (T)                300-1200 K
Pressure (p)                               1-100 (-250) bar
Air-Fuel equivalence ratio ([lambda])      0.6-1.7(-flammability limit)
Residual exhaust gas (EGR)                 0-50 mass-%
Secondary component ([H.sub.2],            0-40 mol-%

Table 2. Adapted calibration parameters for pure methane

Parameter          Value                      Unit

[Z*]             0.0550                -
[E.sub.t]           57961                     K
[B.sub.i]               1.22878*[10.sup.18]   bar
m                       1.5                   -
r                       0.985                 -
n                       2.439                 -
F                       0.2759                cm/s
G                  -11428                     K
[n.sub.EGR]             1.1934                -
[n.sub.a]               0.8809                -
C                       0.8451                -

Table 3. Adapted spline values

Z*        [lambda]   [S.sub.1](Z*)   [S.sub.2](Z*)

0.08849     0.6        0                      1.01982
0.07682     0.7        0                      1.01286
0.06787     0.8        0                      1.00288
0.06079     0.9        0                      1.00004
0.05504     1.0        0                      1.00000
0.05029     1.1        0                      1.00173
0.04630     1.2        0                      1.00116
0.04289     1.3        0.01807                0.99045
0.03995     1.4        0.02959                0.98325
0.03738     1.5        0.04197                0.97519
0.03513     1.6        0.05563                0.96593
0.03313     1.7        0.07079                0.95531
0.03135     1.8        0.08397                0.94653
0.02975     1.9        0.09715                0.93774
0.02830     2.0        0.11033                0.92896

Z*         [S.sub.3](Z*) [K]     [S.sub.4](Z*)

0.08849     1681.2                     0.6479
0.07682     1788.5                     0.6897
0.06787     1892.0                     0.7179
0.06079     2010.6                     0.6913
0.05504     2117.0                     0.5837
0.05029     1968.1                     0.6635
0.04630     1841.4                     0.7107
0.04289     1732.0                     0.7407
0.03995     1636.6                     0.7616
0.03738     1552.4                     0.7768
0.03513     1477.6                     0.7884
0.03313     1410.7                     0.7978
0.03135     1350.3                     0.80993
0.02975     1295.6                     0.82202
0.02830     1245.8                     0.83411

Table 4. Calibration parameters for alkane expansion

Parameter                                   Ethane    Propane   n-Butane

[m.sub.F]                                    0.4107    0.235     0.4156
[mathematical expression not reproducible]  -0.3409   -0.4075   -0.6320
[mathematical expression not reproducible]   0.2332    0.3434    0.5157
[mathematical expression not reproducible]  -0.4908   -0.5897   -0.7520
[mathematical expression not reproducible]   0.5296    1         0.8018

Table 5. Additional calibration parameters for [H.sub.2]-admixture

Parameter     Value

[a.sub.p]      6.45
[n.sub.p]     -0.10216
[a.sub.lam]   -0.2517
[b.sub.lam]    0.4359
COPYRIGHT 2017 SAE International
No portion of this article can be reproduced without the express written permission from the copyright holder.
Copyright 2017 Gale, Cengage Learning. All rights reserved.

Article Details
Printer friendly Cite/link Email Feedback
Author:Hann, Sebastian; Urban, Lukas; Grill, Michael; Bargende, Michael
Publication:SAE International Journal of Engines
Date:Apr 1, 2017
Previous Article:Development and Validation of a Quasi-Dimensional Dual Fuel (Diesel--Natural Gas) Combustion Model.
Next Article:Phenomenological Autoignition Model for Diesel Sprays Using Reduced Chemical Kinetics and a Characteristic Scalar Dissipation Rate.

Terms of use | Privacy policy | Copyright © 2019 Farlex, Inc. | Feedback | For webmasters