A Physics and Tabulated Chemistry Based Compression Ignition Combustion Model : from Chemistry Limited to Mixing Limited Combustion Modes

Résumé — Un modèle de combustion à allumage par compression basé sur la physique et la chimie tabulée : des modes de combustion contrôlés par la chimie jusqu’aux modes contrôlés par le mélange — Ce papier présente une nouvelle approche 0D phénoménologique pour prédire le déroulement de la combustion dans les moteurs Diesel à injection directe pour toutes les conditions d’utilisation usuelles. Le but de ce travail est de développer une approche physique en vue d’améliorer la prédiction de la pression cylindre et du dégagement d’énergie, avec un nombre minimum d’essais nécessaires à la calibration. Les contributions principales de cette étude sont la modélisation de la phase de pré-mélange de la combustion et une extension du modèle pour les stratégies d’injections multiples. Dans ce modèle, le taux de dégagement d’énergie dû à la combustion pour la phase pré-mélangée est relié à un taux de réaction moyen du carburant. Ce taux de réaction moyen de carburant est évalué à l’aide d’une approche basée sur un taux de réaction local de carburant tabulé et la détermination d’une fonction de densité de probabilité (PDF) de la fraction de mélange (Z). Cette PDF permet de prendre en compte la distribution de richesse existante dans la zone pré-mélangée. L’allure de cette PDF présumée est une β-fonction standardisée. Les fluctuations de la fraction de mélange sont décrites avec une équation de transport pour la variance de Z. La définition standard de la fraction de mélange, établie dans le cas de flammes de diffusion, est ici adaptée à une combustion pré-mélangée de type Diesel pour décrire l’inhomogénéité de la richesse dans le volume de contrôle. La chimie détaillée est décrite au travers de la tabulation du taux de réaction relatif à la flamme principale et du délai d’auto-inflammation relatif à la flamme froide, ces tabulations sont fonction de la variable d’avancement c, du taux de gaz brûlé ainsi que des grandeurs thermodynamiques telles que la température et la pression. Le modèle de combustion de diffusion est basé sur l’évaluation d’une fréquence de mélange qui est fonction du taux de turbulence dans la chambre et fonction de la masse de carburant vapeur disponible. Ce modèle de combustion est l’un des sous-modèles inclus dans le modèle thermodynamique de chambre de combustion basé sur une approche 2 zones différentiant gaz frais et gaz brûlés. De plus, un sous-modèle concernant la multi-injection a été développé pour prendre en compte les interactions entre les jets et pour décrire l’impact de la multi-injection sur le mélange air/carburant. Les résultats numériques ont été comparés aux mesures réalisées sur un moteur Diesel Renault de 2 litres. Pour l’ensemble des points de fonctionnement explorés, une bonne correspondance a été obtenue entre les résultats issus des simulations et les mesures expérimentales. Oil & Gas Science and Technology – Rev. IFP Energies nouvelles, Vol. 66 (2011), No. 5, pp. 823-843 Copyright © 2011, IFP Energies nouvelles DOI: 10.2516/ogst/2011138

phenomenological approach to predict the combustion process in multi injection Diesel engines operated under a large range of running conditions.The aim of this work is to develop a physical approach in order to improve the prediction of in-cylinder pressure and heat release.Main contributions of this study are the modeling of the premixed part of the Diesel combustion with a further extension of the model for multi-injection strategies.In the present model, the rate of heat release due to the combustion for the premixed phase is related to the mean reaction rate of fuel which is evaluated by an approach based on tabulated local reaction rate of fuel and on the determination of the Probability Density Function (PDF) of the mixture fraction (Z), in order to take into consideration the local variations of the fuel-air ratio.The shape of the PDF is presumed as a standardized β-function.Mixture fraction fluctuations are described by using a transport equation for the variance of Z.The standard mixture fraction concept established in the case of diffusion flames is here adapted to premixed combustion to describe inhomogeneity of the fuel-air ratio in the control volume.The detailed chemistry is described using a tabulated database for reaction rates and cool flame ignition delay as a function of the progress variable c.The mixing-controlled combustion model is based on the calculation of a characteristic mixing frequency which is a function of the turbulence density, and on the evolution of the available fuel vapor mass in the control volume.The developed combustion model is one sub-model of a thermodynamic model based on the mathematical formulation of the conventional two-zone approach.In addition, an extended submodel for multi injection is developed to take into account interactions between each spray by describing their impact on the mixture formation.Numerical results from simulations are compared with experimental measurements carried out on a 2 liter Renault Diesel engine and good agreements are found.

INTRODUCTION
During recent years, considerable improvements of High Speed Direct Injection Diesel Engine technology have been achieved, with a strong increase in fuel economy and a remarkable reduction of emissions and combustion noise.These improvements have been achieved with the introduction of complex system-layouts involving a large number of devices such as EGR system or Common Rail Direct Injection System.New strategies in Diesel combustion tend to use HCCI (Homogeneous Charge Compression Ignition) or LTC (Low Temperature Combustion) combustion modes more and more.Indeed, for Diesel applications, the goal is to achieve the combustion of a homogeneous mixture that produces lower soot and NO x emissions than that produced by diffusion flame.Many authors [1][2][3][4] have proven that the combustion of a premixed lean charge can produce very low NO x and smoke emissions.Many strategies are developed to control HCCI combustion such as massive EGR combined with multiple injection events.In [5] authors investigate the potential of high EGR rate with an early pilot injection and a main injection close to TDC (Top Dead Center).These latest developments in engine technology involve strong developments of engine models to describe the effects of fuel injection systems on combustion process.3D detailed engine simulations using CFD can provide detailed results reproducing the in-cylinder combustion process but, due to long execution times, are a rather limited proposition for predictive investigations and parameter studies.Alternatively, the key issue with spatially simplified 0D-models is the degree of accuracy that can be achieved.Improvement of physical representative capability while keeping reasonable CPU performances in order to be embedded in a full engine simulator is a challenging and relevant topic in 0D model development.Phenomenological 0D engine models are of particular relevance to perform a large number of numerical tests at costs that are much lower than those associated with experiments.Recent models distinguish premixed combustion from diffusion combustion.In Barba et al. [6], the premixed part of combustion is described by a turbulent flame propagation formulation, and the diffusion combustion is described as a single stage mixing controlled process.Chmela et al. [7] have proposed a model that relies on the concept of mixing controlled combustion based on the fact that, in today's DI Diesel engines, the rate of heat release is mainly controlled by the instantaneous fuel mass present in the cylinder charge and by the density of turbulent kinetic energy.Chmela and Orthaber [8] distinguish the fuel which burns in the premixed phase from the fuel which burns in the diffusion phase, by using a transition function.The latter models use a great number of empirical coefficients that require experimental measurements to be adjusted, resulting in a lack of predictive capability to accurately describe current challenges such as quasi HCCI combustion, premixed lean charge combustion, multi injection strategies or combustion with high EGR rate.Some phenomenological engine models have been designed to be directly used as a tool for car manufacturers in order to involve engine simulation in engine control design [9] or to explore alternative combustion systems in Direct Injection Diesel engines by providing useful information on spray penetration [10], by computing heat release during combustion with a detailed description of the Diesel jet [11], or by taking into account spray interaction with a swirl [12].To satisfy both the accuracy and time requirements, the Probability Density Function (PDF) concept has been introduced together with some phenomenological spray models to consider the effects of spatial inhomogeneities, such as equivalence ratio distribution, on the combustion characteristics [13,14].For reliable simulation of the initiation of the combustion by Auto Ignition of the fuel and the following pre-mixed combustion, detailed chemical kinetic models have been introduced together with combustion models [15].
The aim of this work is to develop a combined physical and chemistry based approach for combustion modeling in Diesel engines.The major challenge for the combustion model consists in precisely describing the overall engine operating conditions.Dec [16] has proposed a conceptual model of DI Diesel combustion which distinguishes premixed combustion from diffusion combustion to account for split injections (pilot injection(s), main injection(s), post-injection(s), etc.).The premixed part of the combustion process is governed by chemistry, while the turbulent mixing of air and fuel controls the diffusion combustion.The relative importance between these two phases varies greatly depending on engine speed and load and is heavily influenced by the injection process, mainly the quantity of fuel injected during the ignition delay.In DI Diesel combustion, premixed combustion is characterized by an air/fuel ratio distribution inside the spray volume [17], more or less homogeneous according to time after injection and mixing velocity [18].The first Auto Ignition site appears inside the spray, where chemical and thermodynamic conditions are the most favorable [19,20].In this paper, the air/fuel ratio stratification inside the spray region is described by a presumed PDF approach [21].The diffusion combustion model is based on the Barba et al. [6] approach, in which the Heat Release Rate is described with a characteristic mixing frequency.
Today high EGR rates are used together with multi injection strategies to lower NO x and soot emissions and also combustion noise [5,22].In zero dimensional models, high EGR rates are taken into account through the thermodynamical model and the Auto Ignition model.But it is very difficult to reproduce physical phenomena occurring during multiple injections.Indeed, a great part of these phenomena is driven by internal engine geometry, in-cylinder temperature and species distribution and is modified by each previous injection [5,22].In this study, an approach with multi virtual zones is chosen to describe interactions between each spray.
In the first part of this paper, each sub-model describing the physical phenomena occurring in Direct Injection Diesel Engine is presented.Using a limited number of operating points the model calibration procedure is then described.Finally, using experimental results from a 2 liter Renault Diesel engine operated under a wide range of operating conditions and parametric variation, the simulated results are compared with the experimental ones.

Overview
In conventional D.I. Diesel combustion, the combustion process is characterized by a high air/fuel ratio distribution inside the combustion chamber.Indeed, in case of single injection, just after the Start Of Injection, the spray volume grows and a mixture zone is created with vaporized fuel and entrained surrounding gas.The time between Start Of Injection and start of combustion is the ignition delay period.The first site of Auto Ignition (AI) appears inside this premixed zone, where chemical and thermodynamic conditions are the most favorable.The evolutions of thermodynamic conditions, due to piston motion, mixing process and heat release at the first site, produce the multiplication of Auto Ignition sites.In opposition to premixed combustion in Spark Ignition (SI) engine, premixed combustion phase in Compression Ignition (CI) engine is not described by a flame front propagation.Total Heat Release Rate in premixed combustion can be considered as the sum of the Heat Release Rate due to combustion at each Auto Ignition site.This phase of combustion is limited by the mass of air/fuel mixture available, produced during the Auto Ignition delay.This chain reaction leads to a very fast and strong heat release, which gives the characteristic noise of Diesel Engines.Fuel injected after AI burns in diffusion combustion mode which is controlled by the mixing process.In this combustion phase, the fuel burns progressively when the vaporized fuel is mixed with the surrounding entrained gas.
In the latest Diesel engine designs, several numbers of devices, such as EGR system or multi injection strategies with Common Rail Direct Injection System, are used to improve pollutant emission performances.To model multi injection strategies, an extended model of a previous work described in [23] has been developed.
Figure 1 shows the basic principle diagram of the combustion chamber in the case of two injections.For numerical simulations, this principle is applied to multiple injections, with interactions between each spray zone.Arrows represent the link and interactions between each sub-model.
The description of the thermodynamic conditions in the combustion chamber is based on a mathematical formulation of the conventional two-zone approach.This zero-dimensional thermodynamic model assumes that at any time during the combustion process, the cylinder volume is divided into two zones, corresponding to burned and unburned gas regions.In each zone the thermodynamic state is defined by the means of thermodynamic properties and the specific heat of each gas component changes according to the approximated formula from the JANAF thermodynamic properties table.Unburned gas composition corresponds to dry surrounding air (N 2 , O 2 , Ar, CO 2 ) with fuel vapor (n-heptane C 7 H 16 ).N-heptane has been chosen because of its AI characteristics which are close to commercial Diesel ones, especially in terms of Cetane number [24].Burned gas mixture composition results from a complete combustion equation in lean conditions (O 2, N 2 , H 2 O, Ar, CO 2 ).Both zones are assumed to have the same uniform in-cylinder pressure.The unburned mixture and burned mixture zones are each treated as separate open systems, with mass exchange and variable composition, and the governing equations are the mass and energy conservation equations and ideal gas equation of state.Mass flow rate in each zone is deduced from a balance equation corresponding to mass transfer through intake and exhaust valves, mass transfer between the two zones due to combustion and fuel injection rate in the unburned gas zone.Instantaneous heat transfers are modeled with a Woschni model [24,25] for DI Diesel Engine.
This quasi-dimensional thermodynamic model incorporates several sub-models to take into account several physical phenomena (turbulence, vaporization, spray and entrained gas mass flow rate), represented in Figure 1.
As shown in Figure 1, each injection leads to a zone implementing a spray, vaporization and combustion submodels.Flow field inside the combustion chamber is described with one turbulence model which takes into account each injection.
This approach makes use of the model developed in [18,26,27] for HCCI combustion, but in the present work only three species are considered: burned gas, unburned gas and fuel.The description of the combustion process is here also simplified by modeling the two combustion phases by their own sub-model.The main evolution presented in this paper concerns the approach developed for multi injection: contrary to [27], here each spray can interact with all previous sprays.More details are given in the next section.

Spray and Entrained Surrounding Gas Model
The spray model is based on the approach developed by Naber and Siebers [28,29], which is a theoretical approach for the evaluation of the penetration of a non-vaporizing ideal liquid spray (Fig. 2).
In [28], a scaling law of the non-dimensional spray penetration according to the non-dimensional penetration time is given: (1) With the dimensionless time and density: t t t .d is the orifice diameter, v f the fuel injection velocity at the orifice exit, and α i the idealized spray angle associated to the ith spray.This model is valid for stationary conditions.However, in ICE, conditions are transient as density in the combustion chamber varies strongly due to piston motion.In order to take into account ambient conditions variation in the combustion chamber, Jaine [30] and Dronniou [31] suggest deriving the previous penetration model, described with Equation (1), to obtain a non-dimensional velocity: (2) This instantaneous velocity, which can be estimated at every time according to the surrounding conditions, allows rebuilding the non-dimensional penetration of the jet thanks to a simple integration: (3) Basic principle of combustion chamber model for a two-injections strategy.Jaine [30] proposes: In addition, Siebers et al. [28,29] show that all spray angle data can be fitted with the following relationship: (4) ρ f and ρ a are respectively the liquid fuel and entrained ambient air densities.The relationship between the idealized spray angle α i and the measured spray angle θ i is given by [21]: (5) Through injection experiments in a constant volume vessel containing nitrogen at given pressure and temperature, correlations (3) and ( 5) have been validated.In order to validate the model of gaseous spray tip penetration, injections of Normafluid (fluid used to simulate commercial Diesel) in a closed pressurized constant-volume vessel were computed.Figure 3 presents some computed results compared with experiments.
Figure 3 compares the computed spray penetrations with experimental data obtained in a high pressure vessel for low and high ambient-air pressure conditions.Results are in good agreement with experiments.As expected the model is sensitive to parameters such as thermodynamic state of ambient gas and injection pressure.

tan
. tan Finally, the spray volume of zone i is given by the geometrical cone definition formed by this jet: (6) The entrained surrounding gas model is based on mass and momentum conservation laws.Spray velocity is assumed decreasing with the surrounding air entrainment in the spray volume.Air + EGR mixture density is assumed to be very close to the total density in the combustion chamber: ρ a ≈ ρ.Thus the total entrained gas mass flow rate into the spray zone can be written as [23]: (7) with V cyl and V S i being the chamber volume and the ith spray volume respectively, V ⋅ cyl results from the piston kinematics and V ⋅ S i results from the spray sub model.Equation (7) gives the total entrained gas mass flow rate in the current spray.As shown in Figure 4, the entrained gas into the current spray zone "i" comes from the ambient gas zone, and also from previous spray zones Figure 4 gives the basic principle of multi injection approach in the case of three injections.As an example the zone i, associated to the ith injection, has the i-1th injection and i-2th injection as source zones for gas entrainment, in addition to ambient gas zone.
In this study, the fraction of the total gas originating from previous injection zones and ambient gas zone in the total entrained gas is taken proportionally to each mass fraction for each zone.(8) with m ⋅ in,S k →S i the gas mass flow rate from zone k = 1, 2, ..., i-1 to the ith zone, m ⋅ in,S i the total entrained gas mass flow rate into the ith zone, given by Equation ( 7) and: the mass fraction of the zone k.Comparison of predicted and measured temporal spray tip penetration and spreading angle for two ambient conditions (blue: ρ a = 17 kg/m 3 , red: ρ a = 50 kg/m 3 ) and other parameters kept constant.The orifice diameter, pressure drop through the orifice and injected fuel temperature were 128 μm, 1 600 bar and 350 K, respectively.Because of fuel spray impingement on combustion chamber walls, experimental penetration lengths are limited to mm.
Schematic of multi injection approach with three injections: repartition of gas mass flow rate between zones.
early in the cycle, the second injection rate corresponds to a pre-injection just before the main injection (3rd injection rate) and the 4th injection rate corresponds to a post injection occurring very late in the engine cycle.Figure 5b shows the gas mass evolution in each zone associated to each injection rate.With this approach the zone associated to the main injection (red) interacts with the two previous zones, corresponding to pilot and pre-injection zones, but not with the ambient gas zone.Simulation results show that this multi injection approach well represents the fact that pilot and pre injections modify the surrounding conditions for the main injection.

Vaporization Model
The vaporization model is used to predict the vaporized fuel mass flow.It must take into account spray expansion, aerodynamics, gas temperature and composition in the combustion chamber.Three approaches are mainly used [6,7,18]: the overall vaporization model, the model based on the "d 2 " law and the isolated droplet model.The first approach was selected because of the low CPU time required.It consists in a global first order model vaporizing the injected liquid fuel m ⋅ liq,i according to a characteristic time τ vap [7]: (9) where: (10) with, k the turbulent kinetic energy given by ( 14), C vap the vaporization model constant, d the diameter of the injector hole, T a,∞ the characteristic temperature of surrounding air.This temperature is assumed to be the mean combustion chamber temperature, T a,∞ = T cyl .Finally, the available fuel vapor mass in the zone associated to the ith injection is defined as: (11) with, m ⋅ vap, i the rate of vaporized fuel obtained with Equation ( 9), and the fuel rate consumed by combustion in the ith zone.m ⋅ f, comb, i, pre and m ⋅ f, comb, i, diff are deduced respectively from Equations ( 17) and (29).m ⋅ f, in, S i is the fuel vapor mass included in the total entrained gas mass in the ith zone computed with Equation (7). Figure 5c shows the fuel vapor mass evolution in each zone associated to each injection rate, for a strategy with four injections in non-reactive conditions.Without combustion, all vaporized fuel associated to pilot and pre injection is entrained into the main zone.

Two States Turbulence Model: Modified k-K Model
In Internal Combustion Engine, internal aerodynamics plays a major role.The dynamics of turbulent motion can be described as follows.Firstly, kinetic energy associated with the mean flow is transferred to the large scales of fluctuating motion.Then, kinetic energy is isentropically transferred from large to small scales of turbulence, and finally kinetic energy is dissipated by viscous friction.
Chmela and Orthaber [7,8] demonstrate that the influence of the spray kinetic energy, compared to other production terms, is dominating and thus assume that the other terms are negligible.Indeed, the production of turbulence is mainly due to a strong shear at large scale between spray and surrounding gas in the combustion chamber.But newer engines implement technologies such as swirl flaps which produce large swirl numbers.Thus, the model must take into account combustion subjected to large swirl numbers that influence average movement in the combustion chamber.Swirl which is used for pollution control (reduction of particulate emissions) affects the mixing rate.Two control volumes are then considered in the energy balance equations, as described in Figure 6.Firstly, in the combustion chamber volume, swirl creates a flow and the associated kinetic energy.After the Start Of Injection, swirl dissipates energy in the spray volume, resulting in spray deformation, and kinetic energy associated with the average movement is increased.In the second volume, there are two main steps.At first, the kinetic energy associated to the mean movement is transferred to the fluctuating movement and creates turbulent kinetic energy at large scale.Then, this turbulent kinetic energy is transferred to small scales and is dissipated by viscous friction.
For total flow in the combustion chamber, the kinetic energy associated with the average movement is [32, 33]: (12) In the above equation, m is the total mass in the cylinder and K cyl production terms are mainly related to the kinetic energy of intake, exhaust flows and swirl.In Equation (12), swirl term is defined as: Diagram of the energy balance in the combustion chamber.and the cylinder inertia moment is defined as J = ∫∫∫ V ρ r 2 dr dθ dz.J is assumed to be the sum of the inertia moment of the bowl volume, which is constant, with the inertia moment of the rest of the combustion chamber, which is variable: Φ bore is the bore diameter, Φ bowl is the bowl diameter, H bowl and H piston are the bowl height and the piston position respectively.Density ρ is assumed uniform throughout the combustion chamber volume and it is assumed that ω air = ω swirl .Finally, only viscous friction induces energy dissipation and with the definition of the friction viscous torque: with C F the torque of viscous frictions calculated with the Stokes law in case of high velocity and C diss the dissipation model parameter.
In the following equations, the control volume corresponds to the spray volume.The assumption for the energy balance equation in the spray volume is that swirl does not induce mass transfer between the surrounding air and the spray volume, but only dissipates energy.By writing the energy conservation equations, the rate of kinetic energy of average and turbulent movement in the spray are [32, 33]: (13) (14) With the production and dissipation terms, P = c P K/t P , D = c D k/t D , and characteristic times: t P = L I /U f , t D = L I /u'.The integral length scale L I is assumed proportional to a characteristic dimension of the combustion chamber, Finally, kinetic energy due to the injection rate is written as: C dish , n hole , S hole , and ρ f are the discharge coefficient, the number of injector holes, the injector hole surface and the fuel density respectively.

Interactions Premixed/Mixed Controlled Combustion
The distribution of the vaporized fuel mass between premixed and diffusion combustion is done with a constant parameter C int ∈ [0,1].When the ratio burned/total injected fuel mass is higher than the parameter, all the vaporized fuel is burned in the diffusion combustion mode.All fuel vaporized, before the parameter is reached, is added to the premixed zone.This definition of fuel vapor distribution allows taking into account a delay between the Auto Ignition in premixed zone and the beginning of the diffusion combustion.Moreover, Barba et al. [6] observe that at the beginning of the diffusion combustion there still is a chemical delay retarding the combustion.To take into account this delay, Barba et al. [6] propose a simple empirical equation: Equation ( 15) is nothing else than a pre-factor for the diffusion combustion, described with Equation (29), retarding the mixing-controlled burning while the premixed combustion is not finished.The value of the exponent is chosen equal to: n = 3.

Premixed Combustion Model
The accurate description of the premixed combustion is important especially at part load conditions or for multi injection strategies.Indeed, for single-injection tests with lowtorque and strongly diluted conditions, the most part of the fuel vapor burns in premixed mode.The premixed combustion influences the diffusion combustion by consuming a great part of the injected and vaporized fuel which is then no more available for all the subsequent combustion phases.The same phenomenon occurs for pilot combustion in the case of multi injection strategies.The fuel mass burned in pilot combustion must be known for an accurate description of the main combustion phase.One of the main issues in conventional Diesel combustion modeling is that Auto Ignition takes place in stratified-mixture conditions [20,34].In this work, the premixed combustion mode is modeled assuming that the reaction mechanism, which controls premixed combustion, is similar to the one that takes place during ignition delay.Consequently, to model this part of combustion, a detailed chemistry-based Auto Ignition model including low temperature phenomena is used to compute a local reaction rate of fuel.The importance of chemical kinetics in new injection strategies in conventional Diesel engines does not allow the use of classical Diesel Auto Ignition models based on oversimplified representations of the chemistry such as those found in Barba et al. [6].The tabulation method adopted in this work has been originally developed for the 3D Internal Combustion Engine calculations ECFM3Z model [35] and here adapted for the 0D approach.
Moreover in DI Diesel engines, an equivalence ratio gradient exists across the spray volume [36].Due to the non-linear dependency of the reaction rate of fuel to equivalence ratio, the equivalence ratio distribution, which develops in the surrounding gas, must be considered to correctly estimate the rate of heat release.Thus, in the present model, the mean reaction rate of fuel is evaluated by an approach based on the determination of the Probability Density Function of the mixture fraction Z i ∈ [0,1].The mean reaction rate associated to fuel is given by the convolution product between the local reaction rate of fuel ω ⋅ f,i , determined by tabulation, and the Probability Density Function of the mixture in the control volume : Equation ( 16) is obtained according to several assumptions [17,26,27]: temperature, pressure and EGR mass fraction are considered homogeneous in the spray volume, and the progress variable c i is assumed homogeneous during combustion, i.e. c i (Z i ) = c ~i, ∀Z i .Thanks to these assumptions the general form of the joint PDF [18,37] is reduced to a simple single variable PDF.In [37], the impact of this assumption on the combustion process is discussed and it is shown that the Auto Ignition process is not affected.The rate of fuel mass consumed by the premixed combustion is then: (17) m f,pre,i is the mass of the vaporized fuel available for premixed combustion in the ith spray zone.This mass is obtained from Equation (11) and with the interaction premixed/diffusion combustion model described previously.

Tabulated Auto Ignition Model Based on Detailed Chemistry
Owing to its low CPU time cost, a tabulation technique issued from three-dimensional approaches [38][39][40] is chosen in this work.The tabulation strategy used is based on the method developed by Colin et al. [39].The method is based on tabulated reaction rates to capture the early heat release induced by Low Temperature Combustion.Cool flame ignition delay when present and cool flame fuel consumption are also tabulated.The reaction rate, fuel consumption, and cool flame ignition delay tables are built a priori from complex chemistry calculations.The Auto Ignition (AI) model has been first validated by Colin et al. [39] in homogeneous configurations and then applied to 3D simulation of conventional Diesel combustion.This approach, previously developed by Pires da Cruz [38], relies on the a priori construction of a database from complex chemistry simulations of Auto Ignition.The Auto Ignition phenomenon can essentially be seen as the rapid oxidation of a premixed charge of fuel and surrounding air.In ICE, Auto Ignition can include cool flame chemistry, depending on local conditions and the main challenge for modeling AI, is to predict when and where this phenomenon occurs.Simulations in homogeneous reactor [41] show the existence of two delays with these operating conditions.A small fraction of fuel is rapidly released after a first delay and then reactions slow down until a second delay is reached.The present model uses a single tabulated delay and a tabulated database for the reaction rates as a function of the progress variable c i .
In [39], the Auto Ignition database has been generated with the SENKIN code [39,42] and calculations have been performed at constant pressure.The Curran et al. [43] mechanism for n-heptane has been chosen as the reference fuel, due to its ignition delay characteristics close to commercial fuels used for Diesel engines.Database points are supposed to cover the widest range of possible combinations of thermodynamic conditions that can be found in conventional Diesel or HCCI engine.
With this Auto Ignition model, the cool flame ignition delay τ LT is calculated according to [38].τ LT is defined as the instant of maximum temperature versus time derivative, (dT/dt) max .Before the low temperature delay is reached, it is very difficult to predict Auto Ignition with reduced chemical models since low temperature chemical reactions are very complex and highly non-linear.The cool flame delay is calculated through integration of the intermediate species Y I .Considering a homogeneous flow without convection and diffusion, the growth rate of Y I is proportional to the amount of fuel tracer Y Tf,i and is a function of the local cool flame ignition delay defined earlier: (18) with: (19) B is a characteristic time chosen equal to 1 s.[38].The fuel tracer represents the amount of fresh mixed fuel that would exist in the spray volume without chemical reactions.The cool flame ignition delay is reached when the intermediate variable is equal to the amount of fuel tracer present in the fresh gases: Y ~I ≥ Y ~Tf,i .Once τ LT is reached, a certain amount of heat c 1 is released.This parameter is tabulated using complex chemistry calculations to obtain a realistic early heat release in all cases.If Y ~f > (1c 1 ) Y ~Tf , the cool flame fuel consumption is computed from: (20) In the original model [38], the fuel consumption characteristic time τ c is assumed constant, but in the present model it is defined as a model parameter.Indeed, this parameter has an impact on the development of the second delay.It is shown thereafter that this parameter can be kept constant for the overall operating range for one particular engine.In the absence of cool flame, both delays are equal and c 1 = 0, so the reaction rate is computed according to the following section.The remainder of the Auto Ignition process is computed , from tabulated results of the reaction rate as a function of the progress variable c i , again issued from complex chemistry calculations.The instantaneous reaction rate has been tabulated a priori for chosen values of c ~i as a function of (T 0,i , P 0,i , φ, X RES , c ~i).In [39], pressure, equivalence ratio and fraction of dilution gases are independent of the fact that combustion occurs or not inside the computational cell.The equation defining the local initial temperature of the mixture, T 0,i = T 0,i (Z), as a function of the mixture fraction Z is: (21) where Ξ -1 is the inversion function of the JANAF thermodynamic table and h is the specific enthalpy of the mixture given by: (22) h f (T u ) is the enthalpy of the fuel vapor at unburned gas temperature and h a (T a ) is the enthalpy of ambient gases.This last enthalpy is given by the energy conservation equation for ambient gas: (23) h u (T u ) and h b (T b ) are respectively, the enthalpy of unburned and burned gases at the unburned and burned temperatures computed with the thermodynamic model.Y b,i is the burned mass fraction in the ith premixed zone.With this definition of the local initial temperature, the temperature stratification due to the air/fuel ratio distribution is considered.The time derivative of the progress variable is stored for discrete values of the progress variable.An adequate refinement is needed in the low temperature region.This corresponds to small values of c ~i, and sufficient accuracy was obtained with seven points [39].The definition of the progress variable must satisfy the following properties: c i must be monotonic with the progress of combustion, c i varies from 0 at the beginning of combustion (when all the fuel is in fresh gases) to 1 when all fuel is in the burned gases.c i must be representative of all main reactions.The main assumption of this tabulation method is that the progress variable representative of the fuel consumption in the engine code is equivalent to a progress variable based on temperature in the complex chemistry code.This assumption, valid for an adiabatic incompressible flow with unit Lewis number, is acceptable in the framework of this study where the relevant information is the fuel consumption leading to a temperature increase.The reaction progress variable c ~i is defined as: , pre, i is the fuel mass consumed by premixed combustion, and m f, pre, total, i is the total mass of vaporized fuel available for premixed combustion.After τ LT , while Y ~f,i < (1c 1 ) Y ~Tf,i , the fuel consumption between the cool flame and the main ignition is a function of the progress variable reaction rate ω ⋅ c and is calculated according to: The local reaction rate of fuel ω ⋅ f,i is deduced from Equation (20) for the low temperature delay and from Equation (25) for high temperature delay. (26) In [39], the model behavior was first tested in different homogeneous configurations, academic configurations at constant volume and, in order to reproduce engine conditions, with a volume variation corresponding to an engine rotation speed.The initial composition was varied to test different thermodynamic conditions corresponding to classical Diesel and HCCI engine operating points.In [39], an excellent agreement is observed with the Curran et al. [43] mechanism for several initial thermodynamic conditions.

Presumed Probability Density Function:
Mixture Model Due to the air/fuel ratio distribution in the spray, the premixed part of combustion in DI Diesel engine is not instantaneous.There exist multi-ignition spots due to local thermodynamic properties evolution.With conventional approaches in 0D modeling, this zone is defined homogeneous [30].With new combustion modes, this consideration is not sufficient to estimate an Auto Ignition delay or a reaction rate.This work makes use of a statistical approach which consists in presuming the shape of the stratified zone with a Probability Density Function.Non homogeneous mixing between fuel vapors and surrounding gas is presumed with a standardized β-function PDF shape.Three characteristic shapes of this function can be distinguished according to the mixing process.Firstly, when fuel vapor and the surrounding air are perfectly unmixed ( Zi "2 = inf), the function corresponds to two peaks at Z = 0 (pure oxidizer) and Z = 1 (pure fuel).Secondly, when an air/fuel ratio distribution exists inside the spray, the function is a continuous curve between Z = 0 and Z = 1.Finally, when air and fuel are perfectly mixed ( Zi "2 = 0), i.e. homogeneous mixture, the function is only one peak centered at the mean value of Z.The mixture fraction definition has been used to describe the equivalence ratio distribution in the spray zone.The mean value of the mixture fraction is defined as a conserved scalar ω ω ω i are respectively the mass fraction of fuel and the mass fraction of oxygen in the ith spray volume, and Y O 2 ,a,nm is the mass fraction of oxygen in the not-mixed zone (ambient gas zone).Zi "2 quantifies the mixture fraction fluctuations, in the gas mixture, related to mixing.Based on the approach used in a previous study [23] and developed by Mauviot [18,37] and Dulbecco [26,27,44], we relate statistical variables to physical phenomena.Deriving the probabilistic definition of the mixture fraction variance, Equation ( 28) is obtained.In this equation we distinguish five contributions influencing the mixture.The first term corresponds to the vaporized fuel flow contribution, the second term corresponds to the entrained gas mass flow from ambient gas zone.Both terms increase the variance.The third term corresponds to the dissipation, it decreases the variance.This term is the micro mixing process contribution and indicates that the mixture state is directly related to turbulence inside the chamber.Within this term, the dissipation rate is related to the integral length scale via, ε = u '2 /L I .The fourth term is related to the entrained gas flow from previous spray zones, this term can increase or decrease the variance according to the mixture state in previous spray zones.Finally, the last term is related to the outgoing flow.This term increases the variance: see Equation (28).
In this equation, m S i is the sum of the gaseous fuel mass plus the mass of entrained air + EGR in the ith spray, m ⋅ vap,i is the vaporized fuel flow given by Equation ( 9), m ⋅ in,A→S i is the entrained surrounding air flow in the ith spray volume from ambient gas zone, given by Equation (8).C Φ i is a constant parameter linked to turbulence that has to be adjusted.m ⋅ in,S k →S i and m ⋅ out,S k →S n are respectively the entrained mass flow from previous spray zone and the outgoing mass flow to the next spray zones, given by Equation (8). Figure 7 gives an example of the mixture state evolution in each zone, corresponding to the case with four injections in non-reactive conditions.
In Figure 7, dashed lines are results without interactions based on the first multi injection approach developed in [23].Figure 7b shows that in the case with interactions, the mean value is higher than without interactions.Indeed, with interactions, vaporized fuel from previous spray zones is entrained in the current zone.Moreover, Figure 7c shows that in the case with interactions the variance is smaller, and then interaction enhances the mixture.In conclusion, this approach well reproduces the main goal of multi injection that is to prepare the mixture state for subsequent injections.

Diffusion Combustion Model
The second part of the combustion model is the diffusion or mixing-controlled combustion model.This combustion mode is based on the approach developed by Barba et al. [6].For engine operating conditions corresponding to full load and with several pilot injections, the main combustion is principally a mixing-controlled flame.This phenomenon, which is controlled by the mixing between surrounding air and fuel vapor, is slower than the premixed combustion.Common zero-dimensional diffusion models are based on the available fuel vapor mass enclosed in the control volume and on a characteristic mixing frequency which is a function of the turbulence density FQ(k).The available fuel mass is the difference between the injected and vaporized fuel mass and the already burned fuel, as defined in the vaporization model.A physical model, which takes into account the rate of injection and the vaporization characteristic time, is then defined.The characteristic mixing frequency, FQ(k), depends on a ratio between the mixing velocity and a mixing length.The mixing velocity is determined starting from the turbulent kinetic energy density in the cylinder, and Equation (14) gives this density evolution.Several solutions were developed to link the mixing frequency with the ratio between the mixing velocity and a mixing length [7,45,46].The approach proposed by Barba et al. [6] is here adopted.
The fuel mass consumed by the diffusion combustion is: φ n N is the number of nozzle holes, φ is the actual fuel/air equivalence ratio, V MP is the mean piston velocity, m f,i,diff is the mass of the vaporized fuel available for diffusion combustion obtained with (11) and according to the interaction premixed/ diffusion combustion model.The two parameters C V MP and C k are optimized using a limited number of characteristic operating points by matching the simulated combustion results with experimental heat release profiles.The determined values are then kept constant for all engine operating conditions.

MODEL CALIBRATION
To evaluate the overall combustion model, operating points covering the whole engine operating range have to be used (Fig. 8).Indeed, the model is evaluated on several steady state operating points covering a large range of engine speeds, loads, EGR rate, and injection strategies representative of the whole engine application field.The various experimental measurements were carried out on a 2 liter Renault Diesel engine with a 16.2:1 compression ratio, equipped with a common rail system.Injection rates are modeled with Hermite polynomials optimized for each injection.
The overall combustion chamber model is calibrated with a limited number of operating points, using the experimental data provided by the engine bench.The determined model parameters are then kept constant for the whole range of engine operating parameters.The global methodology consists in separating the different sub-models calibrations.Therefore in order to limit the interactions, each model is calibrated on specific operating points.To test the model predictivity a minimum set of engine tests is chosen, less than 10% of all available tests.Turbulence sub model parameters were identified on 3D computation results.These parameters in Equations (12,13) and ( 14) are chosen constant for all simulations.Spray sub model parameters were identified on experimental results performed with a Bosh injector in a closed pressurized constant-volume vessel.Same injector type was used in the 2 liter Renault Diesel engine.Some results are shown in Figure 3.
Theoretically, each injection is characterized by six parameters.These parameters are summarized in Table 1.
It is assume that in pilot and pre injection all vaporized fuel burns in premixed mode.Thus, only three parameters for pilot and pre injection are necessary: vaporization sub model parameter and premixed combustion sub model parameters.Only twelve parameters have to be adjusted for a strategy with three injections.All subsequent results are presented, with this set of parameters kept constant.

Analysis of an Operating Point
In this section the operating point, described in Table 2, is studied in detail.
Figure 9a shows experimental and computed in-cylinder pressure.As can be seen, the computed results are in good agreement with experiments.Moreover, computed in-cylinder temperature is also in good agreement with temperature computed from experimental pressure (Fig. 9b).Both temperatures are evaluated by taking into account a variable mixture composition given by the model.Simulated in-cylinder temperature is computed from unburned and burned gas temperatures.Thus it can be assumed that gas temperatures are well estimated and could be used for pollutant formation models, especially for NO x and CO modeling.
In Figure 10, dashed lines are results without interactions between sprays, based on the first multi injection approach developed in [23].Figure 10a shows the variance evolution in each zone with interaction and without interaction.This figure shows that variance for main injection is lower with interaction than without interaction which means that spray interactions enhance the mixture formation.
Figure 10b shows that the mean mixture fraction, corresponding to the pilot zone, is constant even after the beginning of the main injection.Indeed, the pilot mean mixture fraction is not influenced by the main injection because the combustion has already occurred in the pilot spray zone (Fig. 10).It also can be seen that the mean value and the variance of Z for the pilot injection are not influenced by the interaction model between sprays, as the curves corresponding to the cases with or without interaction are mainly superimposed.For this operating point, when the main injection occurs, the pilot zone is characterized by an almost homogeneous mixture, representative of a low value of the variance of Z (Fig. 10a).In the last term of Equation (28), the outgoing mass flow is computed using the mean value of the mixture fraction, which is representative in this case of the almost mixture state of the pilot zone.Moreover, the variance of Z being low, the term related to the outgoing flow in Equation (28) has a negligible influence on mixture fraction fluctuations.Consequently, the outgoing mass flow from the pilot zone to the main zone doesn't modify the mean value of Z neither the variance for the pilot zone.
The maximum of the pilot mean mixture fraction is equal to 0.01, which corresponds to an equivalence ratio of 0.2.Thus, the combustion for pilot injection occurs in very lean conditions.In this multi injection strategy, the main injection interacts mainly with the pilot zone.Composition of the entrained gas in the main spray zone is then close to the pilot zone composition.Moreover, combustion in the pilot zone has occurred, therefore the proportion of entrained burned gas in the main zone is higher and the proportion of entrained oxygen in the main zone is smaller than in the case without interaction.Thus, equivalence ratio for the main injection is higher in the case with interaction.
Figure 11 shows the Heat Release Rate (HRR) associated with the combustion process.HRR in Figure 11a is calculated from in-cylinder pressure traces and with the heat capacity ratio given by the model.As can be observed in Figure 11a the computed Heat Release Rate is in good agreement with experimental result.The model tends to overestimate the peak for cool flame in pilot combustion (first peak).This is probably due to the assumption that the combustion progress variable can be described at each time using the mean value of the progress variable, which is considered independent of the mixture fraction: c i (Z i ) = c ~i, ∀Z i .Figure 11b shows HRR for premixed and diffusion combustions associated to pilot and main injections.It can be observed that the vaporized fuel corresponding to the pilot injection burns totally in premixed-flame conditions, while the vaporized fuel associated to the main injection burns mainly in Results for the studied operating point described in Results for the studied operating point described in Table 2: a) variance of mixture fraction in each zone with and without interactions; b) mean value of mixture fraction in each spray zone with and without interactions.Apparent heat release rate (J/°) Results for the studied operating point describe in Table 2: a) computed and experimental Heat Release Rate; b) computed Heat Release Rate for premixed and diffusion combustion with interactions and without interactions (dashed lines).
diffusion-flame conditions.Furthermore, Figure 11b shows that in the case with interactions, HRR for premixed combustion associated to main injection is greater than without interaction.

Overall Operating Range
To evaluate the overall combustion model, operating points covering the whole engine operating range have to be used (Fig. 7).Indeed, the model is evaluated on several steady state operating points covering a large range of engine speeds, loads, EGR rates, and injection strategies representative of the whole engine application field.
The following figures present results for all tests, summarized in Figure 7, comparing the most important parameters in combustion prediction: peak pressure value, peak pressure position, IMEP values and crank angle values at 5%, 10% and 50% of burned fuel.
In Figure 12, tests in the tolerance range are represented in blue and tests out of the tolerance range in red.For all operating tests a good agreement is obtained, indeed for maximum pressure 90% of tests are in the tolerance range and for all other combustion parameters more than 70% of tests are in good agreement.It is important to notice that the same set of parameters has been used in the combustion model for all of these tests and that less than 10% of available tests have been used to tune the parameters.The difference between the computed and measured values of IMEP is mainly due to a not Predicted vs measured combustion parameters for all tests resumed in Figure 7, dashed black lines represent the maximum allowed values and black lines represent the perfect compatibility between predicted and measured parameters values.Tolerance ranges are, for Pmax ±2 bar, for alpha Pmax ±1°CA, for IMEP ±1 bar, and for CAx ±2°CA.
accurate modeling of the wall losses during expansion at the end of the engine cycle.
CA05 figure shows that some points are out of the tolerance range.For some of these tests, the crank angle value when five percent of fuel is burned is overestimated.This bad agreement is due to the definition of the combustion progress variable that is described at each time using a homogeneous progress variable based on the amount of fuel consumed, which can lead to delay the first stages of combustion by averaging this input parameter of the reaction rate.
For some of the tests out of tolerance, the crank angle value, when five percent of fuel is burned, is underestimated.These tests correspond to late main injection timing, resulting in a possible interaction between spray and wall which is not taken into account by the model.Due to this wall interaction, the model has a bad evaluation of the vaporization rate and entrained surrounding gas rate, the available mass for combustion is then overestimated.
CA50 figure shows that some points are out of the tolerance range and for the majority of these points the crank angle when 50% of fuel is burn is overestimated.These results show that the end of combustion is not accurately modeled, mainly due to wall interactions, which are not taken into account in the presented approach.Combustion parameters for multi-injection tests with EGR rate variations and with SOI variations (operating conditions: 1 800 rpm, MEP = 5 bar, 11% < EGR < 30% and -6°aTDC < SOI < +2°aTDC).Dashed black lines represent the maximum allowed values and black lines represent the perfect compatibility between predicted and measured parameters values.Tolerance ranges are, for Pmax ±4 bar, for alpha Pmax ±2°CA, for IMEP ±1 bar, and for CAx ±2°CA.are used.As previously mentioned, these results have been obtained with one set of parameters.In spite of a difference on the values of CA50, the model results well follow the trends.The difference between the computed and measured values of CA50 is due to a bad modeling of the end of diffusion combustion, mainly due to the increase of burned gases in the combustion chamber and wall interactions.These results show that the model is adapted to engine pre-mapping as it is able to describe the combustion process with parametric variations.

CONCLUSIONS
This paper presents an original approach for Diesel engine combustion quasi dimensional modeling.Main contributions of this work are the new description of the premixed combustion phase and the extension of the model to multi injection cases.The premixed combustion model makes use of a statistical description of the air/fuel ratio distribution inside the spray region.The extended model for multi injection is based on a repartition of entrained gas flow into the current spray zone.Moreover, this extended model takes into account the impact of mixture state in previous spray zones on the mixture state in the current spray zone.
Operating points with large mixing times have been used for model validation and a good agreement between experimental and computed values has been shown.For a wide range of engine operating points, simulated cylinder pressure and Heat Release Rate traces show a good agreement with experimental measurements.The overall objective being the ability of the model to reproduce combustion modes ranging from conventional Diesel mode to full HCCI mode, future work should improve the model to take into account in-cylinder temperature distribution to model full HCCI combustion.To improve the existing model a more physical vaporization model is required, especially to model cold start conditions.Indeed, in this case, vaporization affects the ignition delay, the shape of the PDF and the fuel mass available for each combustion mode.Furthermore, a more detailed spray model can improve general results especially in cases with multi injection.Indeed, wall impingement has a great impact on mixing process and combustion.

Figure 5 N
Figure 5 represents a general case with four injections in non-reactive conditions.As shown in Figure 5a, each injection rate has a typical profile corresponding to the injection type.The first injection rate corresponds to a pilot injection Y m m k

Figure 5
Figure 5Model results for 4 injections in non-reactive conditions: a) injection rate profile; b) gas mass evolution in each zone; c) fuel vapor mass evolution in each zone.

Figure 7
Figure 7Model results for an example with four injections in non-reactive conditions: a) injection rate profile; b) mean value of mixture fraction according to crank angle degree in each zone; c) variance of mixture fraction according to crank angle degree in each zone.

Figure 8
Figure 8 Set of experimental engine tests representative of the overall engine operating range (347 tests).In this figure, engine tests are identified with engine speed and IMEP and also with EGR rate and injection strategy.

3. 2 . 2
Figures 13 shows results for pressure peak value, pressure peak position, IMEP values and crank angle values at which 5%, 10% and 50% of fuel is burned.These results are obtained for 12 multi-injection tests with EGR rate variations and with SOI variations.For each SOI, 3 values of EGR rate

Physics and Tabulated Chemistry Based Compression Ignition Combustion Model: from Chemistry Limited to Mixing Limited Combustion Modes
-This paper presents a new 0D et al. / A Physics and Tabulated Chemistry Based Compression Ignition Combustion Model: from Chemistry Limited to Mixing Limited Combustion Modes

Table 2
Bordet et al. / A Physics and Tabulated Chemistry Based Compression Ignition Combustion Model: from Chemistry Limited to Mixing Limited Combustion Modes : a) computed and experimental in-cylinder pressure; b) computed and experimental in-cylinder temperatures, burned and unburned gas temperatures from simulations.N