Development of a FPI Detailed Chemistry Tabulation Methodology for Internal Combustion Engines

Résumé — Développement d’une tabulation FPI pour la prise en compte de la chimie complexe dans la simulation 3D moteurs — L’objectif de cette étude est d’utiliser une tabulation de la chimie de type FPI ( Flame Prolongation of ILDM ) pour la simulation 3D de la combustion dans les moteurs. La première di ﬃ culté a été d’adapter la méthode à des codes compressibles tout en se limitant à des tailles de tables raisonnables. Pour satisfaire ces contraintes, une nouvelle formulation a été proposée. Elle permet de garder une structure de code basée sur le transport de quelques espèces chimiques tout en assurant la cohérence de l’évolution du système grâce à une équation pour l’avancement de la réaction. Le terme source chimique impliqué dans cette dernière est directement extrait de la tabulation. L’approche retenue permet également d’appliquer la modélisation FPI à des problématiques moteur comme l’utilisation de carburants complexes caractérisés par des chimies détaillées très lourdes ou la combustion fortement diluée. Le modèle proposé, Abstract — Development of a FPI Detailed Chemistry Tabulation Methodology for Internal Combustion Engines — In this paper, the FPI (Flame Prolongation of ILDM) approach for chemistry tabulation is applied to 3D internal combustion engine simulations. The ﬁrst issue is to adapt the method to a fully compressible solver with a limited database size. Tabulation reduction could be a solution for reducing the database size but, for the present applications, the compositions together with local thermodynamic conditions are varying over very large ranges. Instead, we propose a novel formulation of the model based on the mass fraction balance equation of main species and on the tabulated chemical source term of the reaction progress variable. Then, the model is extended to engine problematics such as dilutant addition and complex fuels handling. The model is integrated in the IFP-C3D internal combustion engine solver and it is validated against complex chemistry calculations in constant and variable volume conﬁgurations. The model is ﬁnally successfully applied to the computation of Diesel engine operating point with spray injection.


INTRODUCTION
Ecological and economical issues motivate engine research to develop advanced technologies and dedicated engine control management.In this context, control targets are becoming more and more difficult to reach.For this purpose, CFD (Computational Fluid Dynamics) is used to better understand new concepts such as Diesel HCCI (Homogeneous Charge Compression Ignition) and LTC (Low Temperature Combustion) combustion modes or Controlled Auto-Ignition (CAI T M ).It is also useful to optimize the engine lay-outs in order to fully benefit from these new combustion concepts.
In particular, one key-point is to reduce pollutant emissions revealing a crucial challenge for combustion and chemistry modeling.Basically, modeling pollutant emissions involve two coupled models: the chemical mechanism and the turbulent combustion model.The chemical scheme needs to be complete enough to recover the correct temperature, pressure and dilution dependence.This chemical mechanism cannot be directly used in the CFD code at least for two reasons: first, it involves untractable computation times; secondly, pollutant chemistry in engines is not usually taking place in a perfectly mixed mixture.Therefore, a turbulent combustion model is necessary to account for the local species and temperature stratification.
The question is then how can we include the complex chemistry information into the turbulent combustion model in a way that maintains reasonable computation times.
Over the past decade, various strategies such as ILDM (Intrinsic Low-Dimensional Manifolds) [1][2][3] or ISAT (in situ adaptive tabulation) [4] have been investigated to answer this question.ILDM [1][2][3] uses chemical characteristic time scale analysis to construct low-dimensional surfaces (manifolds) in phase space for given thermodynamic conditions (pressure, fuel/air equivalence ratio, ...) considering eigenvalues and eigenvectors of the chemical system.This approach allows to parametrize the evolution of the whole chemical system using a small number of coordinates (relevant species mass fraction) that are stored in a look-up table.The details of this procedure can be found in [1,2].The ISAT approach [4] was originally developed for efficient computation of combustion chemistry.It corresponds to a storage/retrieval technique for the calculation of composition changes due to chemical reaction in order to optimize chemical CPU time.The tabulation is built on the fly within the 3D CFD code because its pre-computation for all possible states would lead to an extremely large and unmanageable database as the simulation proceeds.When a constitutive response to an integration point is needed, the system either returns a tabulated response from the database or runs the underlying complex constitutive model adding the response to the database.The accuracy of the method depends on a parameter called error tolerance which controls the retrieving error from the ISAT constructed tabulation.
Unfortunately, the importance of chemical kinetics on engine combustion description does not allow the use of reduced mechanisms such as the Shell model [5], increasing the difficulty to adapt ILDM and ISAT approaches to 3D-engine calculation.Indeed, in order to correctly represent auto-ignition, flame propogation and extinction but also complex pollutant formation, more sophisticated detailed chemical kinetics must be considered involving several hundreds of species and elementary reactions [6,7].Consequently, highly reduced techniques are needed to target manageable tabulation size.To provide such tabulations, two fully similar methods based on flamelet calculations, called FPI (Flame Prolongation of ILDM) [8][9][10][11][12] and FGM (Flame-Generated Manifold) [13] were recently proposed.Originally, the FPI method was based on the tabulation of a set of unstrained premixed laminar flames but it has been extended to diffusion flames [14,15] as well as auto-ignition [16] phenomena.In all cases, the strategy is to represent the chemical path (for given thermodynamic conditions) through a unique parameter: the progress variable which goes from zero in the unburnt gases to one in the fully burnt gases.Hence, combustion characteristics are uniquely related to this progress variable leading to a highly reduced database.
Even if different examples have demonstrated the ability to apply the FPI approach to laboratory and industrial situations [9,14,16,17], its implementation within a 3D engine simulation code has not yet been performed to our knowledge.Indeed, the large number of variable thermodynamic parameters (temperature, pressure, fuel/air equivalence ratio, dilution, ...) in addition to their large range of variation makes the use of reasonable database size difficult to achieve.
A new formulation of the FPI method is proposed in this paper allowing application of this strategy to internal combustion engine industrial calculations.First, we point out the coupling with a fully compressible solver.Then, atomic equations are written to ensure mass conservation of the method for complex engine fuel chemistry.Finally, the method is applied to real engine calculations (including huge ranges of variable operating conditions during the whole cycle) leading to an adapted tabulation generation and methodology applicable to internal combustion engines.
The aim of this paper is to present adaptations and validations of the FPI tabulation method for engine applications which allow to use this approach with a maximum coherence and minimum requirement on the database size.In particular, a serious difficulty is to limit the number of discretization points for the progress variable.
This paper is organized as follows.First, an approach is proposed to use the FPI model within a fully compressible solver which is often the case of internal combustion engine codes.It includes a strategy to deal with complex engine fuel in particular atomic mass conservation issue thanks to specific species reconstruction.The second section is devoted to the definition of a progress variable compatible with variable dilution fractions.Then, the FPI tabulation generation as well as an adapted methodology for internal combustion engine are detailed.Finally, the proposed model is implemented within the turbulent combustion model ECFM3Z (3-Zones Extended Coherent Flame Model) [18] in the IFP-C3D code [19].IFP-C3D is a 3D CFD code using unstructured meshes and ALE (Arbitrary Lagrangian Eulerian) finite volume method dedicated to internal combustion engine simulations.The FPI approach is first validated in constant volume homogeneous reactor configurations for different conditions.It is then extended to variable volume in order to test tabulation crossing interpolation.Finally, it is applied to a Diesel engine case to illustrate its ability to reproduce a real engine configuration.

Relaxation Method Towards Tabulated Manifolds: Database Size Reduction
The FPI method is based on manifold trajectories in the continuity of ILDM, and was initially developped from simple structure flame collection [9][10][11]14].It assumes that the species trajectories (for given thermodynamic conditions) are uniquely related to a progress variable, c, allowing for chemistry representation with a small database.Hence, for each property: where A is the chemical state (species composition, temperature...) of the considered system and − → X is the phase space vector of thermodynamic conditions including, for instance, local temperature, pressure, fuel/air equivalence ratio, dilution, heat losses, ... Two different strategies can be used to couple this tabulation approach with a fully compressible flow solver: 1.The mean density equation is deduced from a balance equation and a Y c (linear combinaison of some species mass fractions) equation is introduced.The Y c combinaison allows the definition of the progress variable c (see Eq. 11).Consequently, species concentrations are directly retrieved from the tabulation and do not need to be transported.The temperature is deducted from an energy equation.In addition, the only reaction rate needed is that of Y c that is also tabulated.2. In other codes, species transport equations need to be solved and temperature is reconstructed from species mass fractions and energy.This is the case for most 3D engine codes because spray evaporation coupling, pollutant modeling and variable density accounting are much easily handeled using species transport equations.In this case, Y c is a function of transported species mass fractions (most of the time the choice Y c = Y CO + Y CO 2 is done).When species are transported (second case above), two different approaches can also be considered to define the reaction rate: -2.1 Species reaction rates, ωi , are tabulated from the FPI collection flamelets, -2.2 Species mass (or mole) fractions are tabulated.
Most of previous works have used option 1, Y c and density equation [12,14] or option 2.1, species transport equations and tabulated reaction rates [20].Because of their target applications involving fewer input data within relatively small variations (in comparison with engine operating points), their number of tabulated points were not strictly limited.So, they were able to perform highly refined tabulations.Unfortunately, in engine applications, numerous operating conditions (temperature, pressure, fuel/air equivalence ratio, dilution...) over a large range limit the number of tabulated points to conserve a manageable database size.In this context, species transport (option 2.1), unlike Y c transport (option 1), can lead to trajectory inconsistencies in composition space if no refined tabulation discretisation is used.
Dealing with unsteady flow configurations, directly allocating species mass fraction is not possible so that species source terms must be considered.In that case, coarse FPI tabulation involves reaction rate integration errors that may be cumulated and, in some cases, induce trajectory divergence in composition space.To overcome this problem, a new approach (option 2.2 above) has been proposed in this paper: species source terms, ωi , are reconstructed from tabulated species mass fraction through the first approximation: where Y i denotes species mass fraction and Y FPI i is the FPI tabulated species.t is the current time and τ is a characteristic time scale allowing to relax towards tabulated chemistry.It must be greater than the current time step but very small compared to the characteristic combustion reaction time.τ is then limited by:  2) is then the corresponding mass fraction for a combustion progress c(t + τ).
The progress variable reaction rate ωFPI c was initially directly taken from the complex chemistry calculation.It was observed in this case that very small increments of the progress variable of the order of 1/100 were necessary to construct the FPI table in order to follow correctly the complex chemistry calculation.Therefore, we propose instead to assume a constant progress variable reaction rate ωFPI c between two consecutive values c j−1 and c j of the FPI table: (5) with time t j such as c(t j ) = c j extracted from the complex chemistry simulations.For the first interval ( j = 1), a linear evolution of the progress variable is imposed assuming it corresponds to a small advancement (c 10 −3 ): This formulation allows us to recover correctly the autoignition delay with no extreme refinement in progress variable close to c = 0. Theoretically, at the beginning of combustion, the progress variable is zero and as well as its source term given by the complex chemistry and no reaction can therefore occur.Instead, Equation ( 6) means that the first tabulated progress variable point (c 1 ) is reached after a physical time t 1 .An alternative strategy is that used in the TKI approach (Tabulated Kinetics Ignition) [22] for the auto-ignition delay, where a fictive species equation is used to recover the moment when the auto-ignition delay is reached.
Even if Equation ( 5) involves an approximation of the evolution of the reaction rate with piecewise constant values (order one in time), it was observed that the proposed methodology (Eq.(2-6)) allows the FPI trajectory to reproduce the detailed chemistry manifold with a reduced number of tabulated points for the progress variable (31 values in our case).Consequently, contrary to the other tabulated variables, reaction rate is not interpolated in the direction of the progress variable.Equation ( 5) rather imposes constant value within each progress variable intervals.
The explained methodology (Eq.2-6) allows FPI trajectories to reproduce the detailed chemistry manifold with a reduced number of tabulated points particularly concerning the progress variable.

Mass Conservation
Complex chemistry for hydrocarbons involves hundreds of species.However, it is not possible to tabulate all of them and a selection has to be made.Let N be the number of species involved in the detailed chemistry.The FPI tabulation may take into account only M species (M ≤ N).If no continuity equation is solved, the combustion energy source term is constructed from the M species source terms: where h f i is the formation enthalpy of the species i.Clearly, M << N species are tabulated because real engine fuels involve hundreds of intermediate species (about 400 species for the DCPR-laboratory [6] mechanism).Consequently, the method does not exactly conserve either mass or energy.M major species must be chosen so that Equation ( 7) should be as close as possible to the actual energy source term Selection of main species is straight forward for light hydrocarbons like methane [23].However, it is more difficult for complex fuels like n-heptane because of the cool flame chemistry involving hundreds of radicals [24].In fact, in spite of their small individual contribution, their total contribution may not be negligible compared to Equation ( 7) approximation.At this point, two different strategies can be adopted: -the mass of species not accounted for in the M species tabulated is reported on N 2 .This method is simple but does not conserve the mass of individual atomic elements (C, H and O), that are "transformed" in N 2 .This solution offers the advantage that only the N 2 profile is perturbed compared to complex chemistry.This method cannot be applied to cool flame chemistry for which omitted intermediate hydrocarbons are crucial to represent correctly the global heat release; -atomic conservation equations are written.For hydrocarbon combustion, H, C and O balance equations are considered.Consequently, three species evolutions are reconstructed and only M − 3 species are retrieved from the tabulation.The selection of the three reconstructed species is rather tricky: as they are reconstructed and represent the mass of all other non tabulated species, their mass fractions differ significantly from that given by complex chemistry.As a consequence, their impact on the combustion energy source term (Eq.7) can be artifically important and lead to approximate evolutions of the temperature.This method is physically more acceptable but it leads to differences between detailed chemistry and reconstructed species profiles.Reconstructed species must therefore be carefully selected.
For the reasons mentioned above, considering engine fuel chemistry, the second method has been kept.
After various tests, we have finally chosen the following set of species: CO, CO 2 , H 2 O, the engine fuel C x H y and the radical H are tabulated from the complex chemistry autoignition simulations.Notice that atomic H is considered because of its energy contribution (see Eq. 7) despite of its negligible mass.To conserve H and O atomic elements, respectively O 2 and H 2 atomic budgets are written.Because most of the energy variation during the cool flame is due to intermediate and radical hydrocarbons, C conservation is the most difficult equation to close.It must involve a hydrocarbon species which energy source term is representative of the global contribution of intermediate species not accounted for.Meeting this requirement is quite complex, and after various tests, considering n-heptane as the fuel, the C 7 H 14 was found to give a good match between the approximate source term of Equation ( 7) and the exact one.Once again, this reconstructed C 7 H 14 species does not represent a real species involved in the reference complex chemistry.It must be seen as a "dummy species" which evolution reflects missing tabulated carbon atoms in the system compared to detailed kinetics.
Finally, for C 7 H 16 with the FPI tabulation of the DCPR [6] where Y i is the mass fraction of species i and Y T i is the concentration of the tracer of species i known from the ECFM3Z framework (for details on tracer species and ECFM3Z model see Reference [18]).M i denotes the molar weight of species i.The reader may notice that the N 2 source term is zero if NOx chemistry is not considered and consequently the atomic budget in N simply leads to no correction of the N 2 mass fraction.As a first step in the FPI-engine development, NOx chemistry is not considered in this study.
Once the five tabulated species (CO 2 , CO, H 2 O, C 7 H 16 and H) have been advanced in time using their reaction rates (Eq.2), the three remaining species (O 2 , C 7 H 14 and O 2 ) are obtained applying the above atomic budgets (in addition to N 2 that is conserved).
Following Equations (8)(9)(10), not only the mass is conserved but also the nature of the atoms which is physically more acceptable than an N 2 adjustment.

CHOICE OF COORDINATES: PROGRESS VARIABLE DEFINITION
As explained previously, the FPI approach assumes the combustion state to be uniquely defined by a progress variable c.In the context of simple hydrocarbon fuels like methane or propane, previous FPI works [9,14,25] have suggested to define the progress variable as a linear combination of major species, Y c = Y CO + Y CO 2 .Nevertheless, concerning complex kinetics involving cool flame mechanism, different mass fraction combinations were proposed in order to better represent cool flame auto-ignition.proposed because of the relatively regular oxygen consumption during the combustion reaction.In each case, the progress variable is defined as : where Y 0 c and Y eq c are respectively the initial and the equilibrium mass fraction of the Y c advancement.
All these definitions have been verified to be monotonic in time (see for instance the progress variable evolution in Fig. 9).Indeed, for the CO + CO 2 definition, no turning points appear even if rich mixtures are considered as illustrated in Figure 1.This means that each value of the progress variable corresponds to a unique set of mixture composition and temperature (Fig. 1).For rich conditions, this phenomenon is attributed to chemical equilibrium between CO and CO 2 as observed in [9].Different progress variable definitions have been tested and the "classical FPI definition" This definition assumes that the initial amount of carbon in the fuel is oxidized into final products CO and CO 2 during the combustion process.However, CO and CO 2 can also be found initially inside the combustion reactor within residual burnt gases or exhaust burnt gases: For instance, considering LTC (Low Temperature Combustion), operating points with large dilution rates (about 50% in mass) are common issues.In that case, at the beginning of combustion, the progress variable must be zero in spite of the fact that Y CO + Y CO 2 is not zero.Consequently, CO and CO 2 coming from combustion reaction must be distinguished from residual burnt gases.For this purpose, the total mass fraction of any species i, Y TOT i is decomposed into a "combustion contribution" Y i and a dilution part Y D i which corresponds to the species mass fraction present prior to combustion: The dilution part, Y D i , can be given by a specific transport equation as explained in [18].In practice, for Diesel engine applications, most of the time, Y D i = 0 except for species CO, CO 2 , H 2 O and N 2 which represent the main constituants of burnt gases.
However, the decomposition coming from Equation (13) assumes that the dilution fraction Y D i is not involved in the combustion process.This is, of course, an approximation since there is no chemical difference nor segregation between dilution gases and combustion products.The authors recognize that making no allowance for residual burnt gases composition can be discussed, in particular in the context of several studies on chemical dilution effects on combustion [28].Nevertheless, this hypothesis is used instead of taking into account complete dilution gases composition because in this case, the complete dilution gases composition would appear as input parameters of the database, which would increase the database size by at least two orders of magnitude.In spite of its simplicity, this assumption is reasonable, as shown by Figure 2, since dilution gases are mainly composed of N 2 with a small fraction of combustion products which are weakly reactive compared to the fuel and oxygen mixture.Consequently, to simplify the tabulation generation, the initial dilution gases are assimilated to N 2 and only their total mass fraction Y D defined by: is considered as an input in the FPI table.Consequently, the progress variable is defined by substracting to the total CO and CO 2 mass fractions, their dilution contribution: This definition corresponds to the "classical FPI" definition assuming initial burnt gases are chemically neutral.It has been verified to be a monotonic function of time for all thermodynamic conditions: it strictly increases in time from c = 0 (fresh gases) to c = 1 (burnt gases).Thus, the combustion reaction state is uniquely defined by this progress variable.

FPI Database Generation
Originally, the FPI method assumes that the flamelet considered for tabulation is a laminar premixed flame or a diffusion laminar counter-flow flame.In the present study we instead consider auto-ignition homogeneous constant volume reactors using the Senkin code from the Chemkin package [29].Although this reactor is not properly speaking a flame, we still employ the expression "FPI table'' to designate it because the methodology is exactly the same as with premixed or non-premixed flames.Homogeneous reactor simulations were chosen because, in the context of Diesel engine calculations, the ignition delay must be evaluated with precision.Premixed or diffusion flames cannot lead to a good estimation of this delay for two reasons: 1.When using premixed flames, an equivalent "ignition delay'' can be obtained by a simple integration of the progress variable reaction rate over time.But this delay does not necessarily correspond to the auto-ignition delay because the space diffusion of species and temperature in the premixed flame alters the chemical path compared to a homogenous auto-ignition.

It was observed that reduced chemical mechanisms
do not allow a precise estimation of cool and main auto-ignition delays over the wide range of pressure, temperature, fuel/air equivalence ratio and dilution rates encountered in Diesel piston engines.Consequently, only complex mechanisms with hundred of species incorporating the most complete description of auto-ignition reactions can be used.These mechanisms cannot be today used in premixed or diffusion laminar flames because of the excessive CPU time required to solve these large systems.
Diesel fuel is known to be composed of many hydrocarbons for which complex chemical mechanisms are not available today.In this paper, n-heptane has been chosen, as a first step, as a surrogate fuel to represent Diesel kinetics.The DCPR kinetics [6] is employed for the complex chemistry calculations.This mechanism includes about 400 species and 1500 reactions and has been widely validated [6,24,30].The mechanism used in this paper was generated automatically by the EXGAS-ALKANES software developped at the DCPR laboratory in Nancy.It includes three sub sets: -A C 0 − C 2 reaction base involving species with up to two carbon atoms.-A comprehensive primary mechanism which only considers initial organic compounds and oxygen as reactants.-A lumped secondary mechanism.The molecules produced in the primary mechanism, with the same molecular formula and the same functional groups are lumped into one unique species without distinction between different isomers.For a complete description of the detailed chemical mechanism used in this work, we refer the reader to [31].
In a piston engine, auto-ignition takes place, most of the time, near the top dead center inducing a strong pressure increase within a small and almost constant volume variation.In this situation, perfoming the a priori chemical kinetic simulations in a constant volume configuration seems more justified than doing so at constant pressure.Consequently, we generated our FPI table from constant volume Perfectly Stirred Reactor (PSR) calculations.

Tabulation of Inputs and Outputs
The full variety of compositions and initial conditions is described by four inputs: the initial temperature, T 0 and pressure, p 0 , the fuel/air equivalence ratio, φ and the dilutant volumic fraction (corresponding to the mass fraction Y D given by Eq. 14).As explained in Section 2, we approximate the dilutant composition by uniquely considering N 2 dilution for the FPI tabulation generation (of course, the real exhaust burnt gases composition can be included in the 3D CFD calculation).The progress variable (Eq.15) is added to describe the combustion reaction state.The following ranges are used to construct the complete FPI database: -8 initial pressures between 10 bar and 100 bar, -6 fuel/air equivalence ratios between φ = 0.3 and φ = 3, -4 dilutant rates between 0% and 80% of the total volume, -55 initial temperatures between 550 K and 1500 K, -31 values of the progress variable c (Eq. 15).Therefore, for all these operating points, auto-ignition is calculated in a constant volume reactor using Chemkin [29] with the DCPR detailed kinetics for n-heptane [6].
However, because during the 3D CFD computation the local thermodynamic conditions vary continuously in space (convection, diffusion, spray evaporation, ...) and time (compression, heat release, heat losses, ...), the definition of FPI input coordinates for each computational cell needs to be adressed.Determining cell composition (fuel/air equivalence ratio and dilution) is trivial.However, local database input temperature and pressure are less easy to find.Indeed, the local temperature T and pressure p are known in each computational cell, but, the FPI table is based on the initial temperature T 0 and pressure p 0 of the constant volume PSR calculation.We now explain how these variables can be linked.
Even if species composition is varying in the constant volume PSR calculation, internal energy is conserved during auto-ignition: with e s the sensible energy of the mixture defined by: where the sensible enthalpy is h s = T T re f C p dT with C p the heat capacity at constant pressure of the mixture and T re f a reference temperature for which the standard reference state is usually used (T re f = 298.15K).
The tracer species (cf.ECFM3Z model [18]) allow to define entirely the initial composition of the mixture Y 0 i .Using Equation ( 16), the initial sensible energy can be calculated: e s (T 0 ) = e− M i=1 Δh f i Y 0 i .This initial temperature T 0 is finally deduced from e s (T 0 ) and Y 0 i , and the initial pressure p 0 is given by the equation of state: where M 0 is the molar mass of the mixture, R the perfect gas constant and ρ 0 the initial density.
Regarding FPI outputs, the progress variable reaction rate ωFPI In the 3D CFD calculation, in each cell and at each timestep, local conditions (T 0 , p 0 , φ, Y D and c) are determined to be used as tabulation inputs.Concerning the interpolation, the method uses polynomial matrix coefficient for each output variable.For instance, dealing with reaction rate, only four interpolation directions are considered (temperature, pressure, fuel/air equivalence ratio and dilution) and no interpolation on progress variable direction is needeed as explained previously.Consequently, the reaction rate is FPI tabulation example for the point: T 0 = 700 K, p 0 = 30 bar, φ = 0.7 and no dilution.Reaction rate tabulation constructed from detailed chemistry using Equation ( 5).Between two tabulated points (progress variable discretization), reaction rate is a constant according to Equation (5).
interpolated using 2 4 = 16 coefficients with: in which T 0 and p 0 are the local cell input temperature and pressure obtained thanks to Equations ( 16) and ( 18), φ is the
local fuel/air equivalence ratio and Y D is the local dilution of each cell.The coefficients a 0 to a 15 are constant only in each hyper-cube defined by p i 0 < p 0 < p i+1 0 , T j 0 < T 0 < T j+1 0 etc. (where i and j suffixes are respectively FPI discretization points in the pressure and the temperature directions).Note that these coefficients a 0 to a 15 are function of the progress variable c.

Methodology for Internal Combustion Engine Applications
The above FPI tabulation technique has been implemented in the CFD engine code IFP-C3D [19].The calculation algorithm of the model is the following (see Fig. 5): -The nine species presented above are transported in the CFD code.Five are tabulated (CO, CO ).-In each computational cell, the amont of dilution gases Y D i is known thanks to the tracer species transport equations.These dilution gases are subtracted from the total species mass fractions as in Equation ( 13).
-At each calculated point, the progress variable is evaluated by Equation ( 15).-Local input conditions (T 0 , p 0 , φ, Y D ) for FPI tables are determined.-Progress variable reaction rate, ωc (T 0 , p 0 , φ, Y D , c i ), with c i < c < c i+1 is obtained by interpolation in the FPI tabulation (along (T 0 , p 0 , φ, Y D ) directions) and the evolution of the progress variable at time t+τ is estimated according to Equation ( 4).-Tabulated species mass fractions Y i (T 0 , p 0 , φ, Y D , c(t + τ)) are retrieved by interpolation in the FPI table.The next time step mass fractions are then deduced: -Reconstructed species are built from atomic balance (Eq.8-10).-Dilution gases are added to recover the total mass fractions (Eq.13).

Constant Volume Reactors
The model is first validated on the simple configuration of a homogeneous constant volume reactor.Considering adiabatic conditions, the FPI-engine simulation and the complex chemistry obtained with Chemkin solver are compared in Figures 6 to 8. Different conditions were tested varying the fuel/air equivalence ratio (from φ = 0.3 to φ = 3) and the initial thermodynamic conditions within engine ranges (from atmospheric pressure to 100 bar, low (< 700 K) and high (> 1200 K) temperature, from no dilution to high EGR (Exhaust Gas Recirculation) level (80%).Figure 6 shows an example of the temperature evolution during auto-ignition of a lean C 7 H 16 /air mixture with Figure 6 Comparison of temperature profiles obtained with the FPIengine model on IFP-C3D code (symbols) and with the Chemkin detailed chemistry computations (solid line).Perfect Stirred Reactor at constant volume for a lean case (φ = 0.7) with initial conditions: T 0 = 700 K, p 0 = 30 bar, no dilution.

Figure 7
Comparison of main species profiles obtained with the FPIengine model on IFP-C3D code (symbols) and with the Chemkin detailed chemistry computations (solid line).Perfect Stirred Reactor at constant volume for a lean case (φ = 0.7) with initial conditions: T 0 = 700 K, p 0 = 30 bar, no dilution.the detailed chemistry of DCPR [6] and for the FPI-engine model run in IFP-C3D.In the same way, Figures 7 and  8 illustrate main species evolutions respectively for a lean (φ = 0.7) and a rich case (φ = 1.5).In all tested cases, the FPI modeling shows good agreement compared to detailed chemistry in terms of reactor temperature evolution but also in terms of species dynamics.It can be noticed that both temperature and tabulated species profiles slightly differ from complex chemistry reference due to the small number (31) of progress variable values.As it can be clearly seen Figure 8 Comparison of main species profiles obtained with the FPIengine model on IFP-C3D code (symbols) and with the Chemkin detailed chemistry computations (solid line).Perfect Stirred Reactor at constant volume for a rich case (φ = 1.5) with initial conditions: T 0 = 700 K, p 0 = 30 bar, no dilution.
in Figure 6, 7 and 8, the temperature and tabulated species evolution is piecewise linear due to the piecewise constant progress variable reaction rate (Eq.5).But by construction, the complex chemistry temperature and tabulated species evolution are exactly recovered at each discretized progress variable point c i .
On the contrary, reconstructed species like O 2 can significantly differ from complex chemistry calculation, especially in the cool flame region (see Fig. 7).As explained above, because all intermediate and radical species are not taken into account (for example: CH 4 , CH 3 , OH, ...), reconstructed species represent all "omitted" species.For instance, oxygen atoms can be involved in many products such as oxygenated intermediate hydrocarbons that are not considered in the present table.Instead, missing oxygen atoms are reported by the oxygen balance equation into O 2 .Consequently, the O 2 evolution is different from that of the detailed chemistry.Authors agree that reporting the oxygen mass balance into the O 2 species might not be the optimum solution since O 2 is a major species involved in pollutant formation for instance.Actually, work is planed to realize a coupling between this FPI model and pollutant models of NOx and soot.In this context, this assumption will be reconsidered and a minor species containing the oxygen element will be selected to close the oxygen atomic budget.Nevertheless, as a first step of modeling development, this approach has been chosen for the present paper.
These first validations show an overall good agreement while considering a very limited number (31) of progress variable discretization points.This will allow us to apply Figure 9 Influence of the characteristic relaxation time τ (Eq.2-4) on temperature description.Line: reference detailed chemistry of C 7 H 16 -mechanism of DCPR [6].Dashedlines: FPI simulation with τ = 0.001 ms, 0.1 ms, and 1 ms.Auto-ignition for a constant volume homogeneous reactor at φ = 1 with no dilution.Initial conditions: T 0 = 700 K, p 0 = 30 bar. the FPI tabulated chemistry approach to complex internal combustion engine simulations.Moreover, it can be noticed that even if only nine species are considered, the method allows to correctly reproduce the temperature evolution (see Fig. 6), not only the final value but also the cool flame heat release (reminding that temperature is computed from transported species mass fraction and sensible energy).These first results confirm that C 7 H 14 allows a correct representation of the impact of omitted intermediate species on the heat release.It can be pointed out that most important errors are observed, of course, for rich conditions in which unburnt hydrocarbons are numerous but results are still acceptable considering the very limited database size.

Influence ot the Relaxation Time, τ
The presented FPI methodology involves a characteristic time, τ, within Equations ( 2) to (4).In practice, this parameter can be used to handle numerical stiffness which might be problematic in case of a sudden change in the FPI manifold.This situation can happen for instance when gaseous fuel is added by spray evaporation or when very large spatial gradients of mixture fraction are encountered.
First, the influence of parameter τ is investigated on a simple homogeneous constant volume reactor with initial conditions: p 0 = 30 bar, T 0 = 700 K, φ = 1 without dilution.Progress variable evolutions with different relaxation times for Equations (2)(3)(4) are plotted in Figure 9.For this case, the global chemical time τ c (measured as the time needed for c to go from 0.01 to 0.99 in the detailed chemistry calculation results) is approximately 0.5 ms.When τ << τ c (0.001 ms in Fig. 9), the FPI trajectory follows very precisely the detailed chemistry one.For τ = 0.1 ms 0.2τ c , the exact trajectory is still approximately recovered.Finally for τ = 1 ms 2τ c , only the cool flame period is correctly reproduced while an important smoothing of the main flame reaction rate is observed.This test shows that choosing a very small value for τ will allow in any case, to follow precisely the FPI trajectory but might lead to reaction rate stiffness.On the contrary, increasing τ may lead to important trajectory errors.For this reason, a limitation of τ (Eq. 3) was introduced which guarranties that τ always remains smaller than a characteristic chemical time.The evolution of the effective relaxation time τ (limited by (1 − c)/ ωc ) is given in Figure 10 for τ = 0.001 ms and τ = 1 ms.In the case τ = 0.001 ms, no limitation occurs excepted close to the maximum reaction rate occuring between times 6.5 and 6.6 ms.In the case τ = 1 ms, limitations occur during cool flame (between times 6.2 and 6.3 ms) and during an important part of the main auto-ignition.
As expected, species evolution also depends on parameter τ (through relation ( 2)).Hence, the smaller τ is, the better species variations are captured as it is illustrated in Figure 11 for τ = 0.001 ms.On the contrary, if τ becomes large, important errors can be found during reaction: for instance, the CO evolution corresponding to τ = 1 ms in Figure 11 shows an important under-prediction during the reaction.It must be emphasized, that, whatever the relaxation time and the discretization in progress variable, this method guarantees a convergence towards the equilibrium (state c = 1) for Influence of the characteristic relaxation time τ (Eq.2-4) on CO mass fraction.Line: reference detailed chemistry of C 7 H 16 -mechanism of DCPR [6].Dashed-lines: FPI simulation with τ = 0.001 ms and 1 ms (Circles).Auto-ignition for a constant volume homogeneous reactor at φ = 1 with no dilution.Initial conditions: T 0 = 700 K, p 0 = 30 bar.
tabulated species (see Fig. 11).This feature is specific to this relaxation method and is not found when tabulated species reaction rates (option 2.1 in the introduction) are used.

Homogeneous Engine: Variable Volume
In order to better reproduce engine conditions, a variable volume test case corresponding to an engine rotation speed of 2000 rpm has been computed.These tests allow to check numerical time integration accuracy in particular concerning interpolation in the FPI tabulation across temperature and pressure directions.Once again, direct comparisons with Chemkin computations are performed at the same operating conditions (initial pressure and temperature, fuel/air equivalence ratio and dilution rate) along with engine volume variation law.The present test case corresponds to the following conditions: initial temperature, T 0 = 750 K and initial pressure, p 0 = 20 bar at 30 cad (crank angle degree) before TDC (Top Dead Center).The example is for a stoichiometric case, φ = 1 with no dilution.The engine-like operating system is summarized in Table 1.Results plotted in Figure 12 compare CO and CO 2 mass fraction evolutions predicted by the FPI approach and the reference detailed chemistry.Once again, an excellent agreement is observed.In particular, it must be pointed out that the cool flame behavior is very well predicted as illustrated in Figure 13.

APPLICATION TO A 3D INTERNAL COMBUSTION ENGINE SIMULATION
Finally, the proposed model is applied to an industrial engine test case.For this purpose, the FPI approach is introduced in the turbulent combustion model of IFP-C3D.

The 3D CFD Code
IFP-C3D [19] is a 3D CFD parallel code developed at IFP based on Navier-Stokes equations resolved on unstructured meshes using ALE (Arbitrary Lagrangian Eulerian) finite volume method.The turbulent combustion model is the 3-Zone Extended Coherent Flame Model (ECFM3Z) described in [18].In this model, the local mixture fraction stratification is described thanks to three mixing zones: a pure air zone, a pur gaseous fuel zone and finally a mixed air/fuel zone in which combustion takes place.In this mixed zone, auto-ignition combustion was originally modeled by the TKI model [22].In the present paper, this TKI model is replaced by the proposed FPI approach.Note that as we assume fuel and air to be perfectely mixed in this mixed zone, no statistical averaging over mixture fraction or progress variable of the FPI table needs to be performed like in the FPI-PCM model [14].Therefore, the FPI model can be directly implemented in ECFM3Z.

Engine Applications
To illustrate the FPI tabulation approach on Diesel engine simulation, a full load engine naturally aspiratd configuration at 1250 rpm with a single spray injection has been calculated.This application is a first step towards full FPI complex chemistry tabulation used for industrial engine research.The purpose of this example is to test the complete methodology on a real engine configuration.Nevertheless, at this stage of the modeling development, only qualitative comparisons are considered since no pollutant models for soot and NOx emissions are included.
Calculations are performed on a single-cylinder engine based on a PSA DW10 DI Diesel engine using a 3D wedge mesh with homogeneous initial conditions at inlet valve closure (146 cad before TDC).Table 2 summarizes the geometrical characteristics of the engine and the operating conditions.
Diesel engine pressure: experimental measurement (circles) and calculated pressure profiles (solid line) given by the FPIengine tabulation method.Case-1.
Diesel engine pressure: experimental measurement (circles) and calculated pressure profiles (solid line) given by the FPIengine tabulation method.Case-2.
that our purpose is to demonstrate FPI approach feasibility on a complete engine calculation, the surrogate fuel is supposed to have the same chemical global formulation as C 7 H 16 including its thermodynamics properties.In contrast, liquid properties entering the evaporation rate calculation of the spray droplets, are treated as a single component lumped fuel as proposed in [32] or [33].Two injector positions are simulated: Case-1 operating point corresponds to a restricted nozzle tip protrusion whereas Case-2 injector position is 1.1mm below in the combustion chamber.The 3D computation using the presented FPI-engine model correctly reproduces the experimental pressure curves (see Fig. 14 and 15) for both cases.
Auto-ignition delay, maximum cylinder pressure and maximum pressure timing as well as heat release energy are globally retrieved.Of course, the pressure profiles in diesel engines highly depend on the considered detailed chemistry because the mechanism defines the auto-ignition delay.Indeed, to be able to correctly predict Diesel engine behavior, the auto-ignition timing must be known accurately.To this purpose, it must take into account fine chemistry response and coupling with other physical variations that can be modified by the internal engine geometry or by using different engine control strategies.Concerning the IMEP (Indicated Mean Effective Pressure), FPI-engine model simulations are very close to measurements.For the first case, FPI modeling leads to IMEP = 7.8 bar while experiment gives 7.6 bar.Similarly, for the second case the difference is very small: 6.7 bar for the measurement versus 6.9 bar for FPI prediction.
Species profiles during the calculation, plotted in Figures 16 and 17, show a general good behavior for all considered species (tabulated or reconstructed).Maximum of CO occurs around 10 cad after TDC, then, it is oxidized in CO 2 during the expansion stroke.The "dummy species" C 7 H 14 represents intermediate carbon species.Even if its mass fraction remains very small (about 0.1% of the total mass), it must be included in the list transported species to reproduce correctly the heat release as explained previously.The final computed CO emission is about 0.32 kg/h for the first case, which is comparable to the measurement, 0.26 kg/h.For the second case (Fig. 17), CO final emission is larger than the first one according to experiment observation.These results allow us to confirm the potential of the FPI method.It must be pointed out that, at this stage of the model development, CO predictions differ from experimental measurements mainly because of the assumption of the surrogate fuel and the lack of information concerning the coupling with other pollutants, in particular soot emissions.

CONCLUSION
A new version of the FPI tabulation chemistry is developed in the context of CFD codes transporting species mass fractions.In order to deal with Diesel type engine applications, the FPI approach is modified in three main aspects: -The premixed laminar flame is replaced by homogeneous constant volume calculations to generate the FPI database.These PSR calculations allow the use of large complex chemical mechanisms (with hundreds of species) and a precise estimation of the auto-ignition delay and heat release rate.-The species reaction rate is not given by the complex chemistry tabulation as in the original FPI approach but is instead deduced from tabulated mass fractions and the progress variable reaction rate.-Transported species are divided into tabulated and reconstructed species: tabulated species are chosen as the reactants and main products while reconstructed species need to be chosen with care depending on the fuel and chemical mechanism.These reconstructed species allow to guaranty atomic conservation in the CFD calculation and to reproduce correctly the heat release rate during the combustion.This model has been implemented in the ECFM3Z combustion model [18] in order to replace the auto-ignition model TKI [22].ECFM3Z is itself integrated in the 3D engine simulation code IFP-C3D [19].The FPI-engine model has been validated through comparisons with complex chemistry calculations performed with the Senkin solver of Chemkin for a set of test cases: auto-ignition in a homogeneous reactor with constant and variable volumes over a wide range of thermochemical conditions.Species evolutions as well as reaction heat release are correctly predicted by the FPI-model.
At last, it has been applied to the 3D simulation of a Diesel engine and demonstrated the model capability to reproduce industrial applications.The first results show good agreement with available experimental data.This approach allows a reasonable computation time as well as a manageable FPI database size thanks to the reduced number of progress variable points allowed by the proposed species reaction rate expression.
Future research topics have been identified and will be the subject of coming publications.These include the coupling between the FPI tabulation and other pollutant models such as unburnt hydrocarbons, NOx emission or soot formation.Another stage of the development will consist in applying the proposed model to compounds more representative of Diesel fuel than pure n-heptane.At last, the model achievement will be tested on a wide range of operating engine points (conventional and HCCI Diesel combustion, classic spark ignition, CAI T M engines, ...) in order to verify its capacity to reproduce engine geometric parametric variations or fuel effects.

Figure 1
Figure 1 Evolution of temperature versus Y CO + Y CO 2 for homogeneous constant volume auto-ignition calculations for different values of the fuel/air equivalence ratio.Initial conditions: T 0 = 600 K, p 0 = 20 bar, no dilution.
Equation (5) and species mass fractions Y FPI i (T 0 , p 0 , φ, Y D , c), i = CO, CO 2 , H 2 O, C 7 H 16 , H are directly read from the PSR stored calculations.An example of a FPI tabulation point is displayed in Figures 3 (species) and 4 (progress variable reaction rate).

Figure 10 Effective
Figure 10Effective relaxation time during the calculation: relaxation time τ limited by relation(3).Auto-ignition for a constant volume homogeneous reactor at φ = 1 with no dilution.Initial conditions: T 0 = 700 K, p 0 = 30 bar.

4) Y i FPI
(t + δ t) 2 , H 2 O, C 7 H 16 and H) and four are built from atomic balance equations (C 7 H 14 , O 2 , H 2 and N 2

TABLE 2 Characteristics
Because FPI-engine methodology involves balance atom equations for none-tabulated species, the choice of a representative surrogate fuel is crucial.The definition of a real Diesel fuel detailed kinetics is not yet possible.Reminding Figure 16 Main species mass fraction profiles.Circles: O 2 .Squares: CO 2 .Diamonds: H 2 O. Up triangles: CO.Right triangles: "Dummy species" for carbon balance, C 7 H 14 .Case-1. Figure 17 Main species mass fraction profiles.Circles: O 2 .Squares: CO 2 .Diamonds: H 2 O. Up triangles: CO.Right triangles: "Dummy species" for carbon balance, C 7 H 14 .Case-2.