Session 6
COMBUSTION
Chairman: G. Raithby
PILOT DISTILLATE IGNITION OF PREMIXED AND HIGH PRESSURE INJECTED NATURAL GAS COMBUSTION
I.S. Choi and B. E. Milton
School of Mechanical and Manufacturing Engineering
The University of New South Wales, Sydney, Autralia 2052
Dualfuel engines have existed for many years as a means of using alternative fuels as a distillate replacement. In general, the alternative fuels with the greatest potential are the gases such as propane, (predominant in LPG), methane (the largest constituent of both natural gas and biogas) and the alcohols, methanol and ethanol. Natural gas, being the most likely medium term fuel is that under study here. All these fuels have very low cetane numbers. They can therefore be used easily in spark ignition (SI) engines but not in compression ignition (CI) types. While medium (truck) size SI engines exist, there are problems in developing such engines for large cylinder volumes and many of the advantages of CI engines are therefore lost. In addition, during a fuel replacement phase, dedicated alternative fuel engines can only be used on limited routes where the fuelling infrastructure exists and an engine which can operate on both fuels is a preferred option. The solution is to dualfuel with the alternative fuel being introduced in lean proportions much as in the SI types but with a pilot distillate spray to provide the ignition. One problem with such a system is that the turndown ratio of the distillate injector from full load, diesel only operation to minimum load, pilot operation is well outside the range of conventional injectors. However, new types of electrohydraulic injectors are now being developed such as that described by Yudanov^{1}, which may help to solve this problem.
The natural gas (or any of the vaporised low cetane number fuels) is normally mixed in the inlet manifold by what is essentially a gas carburettor although, to overcome defects of carburetted natural gas systems, some developments are considering direct injection of the gaseous fuel into the cylinder. While the former is the simplest and most convenient method, it has been noted by, for example, Karim^{2} that high precombustion gas levels can lead to endgas autoignition type knock similar to that which exists in SI engines in addition to the diesel knock in CI engines which results from the ignition delay of the pilot fuel. The latter is generally thought to be exacerbated by the addition of the gas (Karim et al.^{3}, Nielsen et al.^{4} and Milton^{5}) although some differences exist on its extent. In addition, the slow burning of methane gas needs consideration in gas injection to vehicle engines and direct gas injection has been restricted to large marine engines.
With all recent development in internal combustion engines, a wide range of numerical models have been of immense importance. The difference between the conventional diesel and the dualfuel engine lies in the mixed combustion process. Models of both the premixed and injected gas dualfuel combustion are therefore required. These developments are reported here. In order to better understand and calibrate these models, constant volume combustion bomb tests have also been carried out. This paper deals predominantly with the models developed for the constant volume process as is now described.
Models used with gasoline (SI) and diesel (CI) engines are either thermodynamic, quasidimensional types as discussed by Heywood^{6}, with coupled phenomenological descriptions of the more important processes or, more rarely, fully multidimensional models. The latter are more demanding computationally but allow greater flexibility in that they can more readily handle such items as the heterogeneous gaseous fuel/air mixture surrounding the diesel spray in a dualfuel engine. Hence, the dualfuel combustion model developed here is based on a simplified multidimensional, finite volume system which provides a compromise between accuracy and computational time.
The core of the process is the diffusion model which describes the time dependent mixing of the gas phase substances (fuel, air and combustion products in differing proportions) throughout the spatial grid which is based here, for simplicity and relevance to the combustion bomb used in the tests, on cylindrical coordinates. The numerical solution was determined using the differential conservation equations of mass, energy and linear momentum, the turbulence of the flow being described by means of the ke eddy diffusivity model of Launder and Spalding^{7}.
REFERENCES
 Yudanov, S.V., Development of the Hydraulically Activated Electronically Controlled Unit Injector for Diesel Engines, SAE Paper 952057, 1995.
 Karim, G.A., Knock in DualFuel Engines, Proc. Inst. Mech. Engrs., Vol. 181, Pt. 1, No. 20, pp. 453466, 1967.
 Karim, G.A., Jones, W. and Raine, R.R., An Examination of the Ignition Delay Period in Dualfuel Engines, SAE Int. Fuel & Lubricants Meeting and Expo., Baltimore, Maryland, SAE Paper 892140, 1989.
 Nielsen, O.B., Qvale, B. and Sorenson, S.C., Ignition Delay in the Dualfuel Engine, SAE 870589, 1987.
 Milton, B.E., Improving the Performance of Small DualFuelled Engines, Paper 22, 1st NGV Int. Conf. and Exhib., Sydney, Australia, 1988.
EXACTFACTORIZATION, COMPACT METHODS FOR
MULTIDIMENSIONAL HEAT TRANSFER: APPLICATION TO COMBUSTION
J. I. Ramos
Departamento de Lenguajes y Ciencias de la Computación
E. T. S. Ingenieros Industriales, Universidad de Málaga
Plaza El Ejido, s/n, 29013Málaga, Spain
It is well known that the use of qmethods for the solution of multidimensional advectiondiffusion
reaction equations such as those occurring in heat and mass transfer, combustion, fluid dynamics,
etc., provides a large system of nonlinearly coupled algebraic equations whose solution requires
iterative techniques. Iterations may, however, be eliminated by employing linearized qmethods which
are based on the discretization of the time derivative and the timelinearization of the nonlinear terms;
however, a large system of linearly coupled algebraic equations still results upon the discretization of
the spatial derivatives. In order to eliminate the coupling between different spatial directions in
linearized qmethods, approximate factorization techniques which reduce a multidimensional problem
to the solution of a sequence of onedimensional ones have been proposed. The linear one
dimensional operators may be easily and inexpensively solved by means of a block tridiagonal solver.
Unfortunately, these approximate factorization techniques introduce a factorization error due to the
factorization of the threedimensional operator into onedimensional ones, and the magnitude of this
error may be very large unless small time and space steps are employed in the calculations.
In this paper, an iterative, exactfactorization method based on linearized qtechniques is proposed.
Since it is an exact factorization method, it does not have approximate factorization errors.
Furthermore, since a linearized qtechnique is employed, the factorized, onedimensional operators in
each spatial direction are linear and may be easily and inexpensively solved. However, the exact
factorization introduces higherorder derivatives and a coupling between the different spatial
directions. This coupling has been accounted for by means of a predictorcorrector strategy as
follows. In the predictor step, an approximate factorization method which introduces factorization
errors is employed, whereas the factorization errors are accounted for in an iterative manner in the
corrector step. Moreover, the mixed fourthorder spatial derivatives have been eliminated by using the
onedimensional operators so that only secondorder derivatives appear in the onedimensional
operators which are solved in an iterative manner. Furthermore, in order to accurately resolve the
steep fronts or boundary layers that may occur, the spatial derivatives have been discretized by means
of threepoint, compact difference expressions which are fourthorder accurate in space, so that the
resulting method is secondorder accurate in time and fourthorder accurate in space, and the
difference equations are formulated in delta form so that, for problems which have steady state
solutions, the asymptotic numerical solution is independent of the time step employed in the
calculations if approximate factorization techniques are employed to obtain the numerical solution.
Since the approximate factorization errors depend on the curvature of the solution, they are expected
to be largest where the curvature of the dependent variables is largest. In addition, since these errors
depend on the time step, the steady state solution, if any, depends on the second power of the time
step.
The use of threepoint, compact difference expressions for the discretization of the spatial derivatives
requires that, in addition to the unknown dependent variables, their first and secondorder spatial
derivatives be considered as unknowns. This implies that, if U is the vector of dependent variables
and n denotes the number of spatial directions, the blocks of the onedimensional finite difference
operators which result upon the factorization of compact, linearized qmethods are 3Nx3N for
advectiondiffusionreaction equations compared with NxN for the case that the spatial derivatives are
discretized with secondorder accurate formula, where N denotes the number of dependent variables.
The large increase in the dimensionality of the matrices is compensated by both the ease and
inexpensiveness with which onedimensional operators can be solved and the values of the first and
secondorder spatial derivatives that they provide. These derivatives may be used to concentrate the
grid points were steep gradients occur. Moreover, the use of the onedimensional operators at three
successive grid points together with the fourthorder accurate finite difference expressions which
relate the first and secondorder derivatives to the discrete dependent variables provides seven
equations for nine unknowns, i.e., the dependent variables and its first and secondorder derivatives;
therefore, the spatial derivatives may be eliminated to obtain a tridiagonal matrix for the dependent
variables.
The exactfactorization, compact method proposed here has been used to study the combustion of a
mixture which is governed by a onestep irreversible chemical reaction in twodimensions, i.e., n =
2, by solving the energy and mass conservation equations, i.e., N = 2, subject to homogeneous
Dirichlet boundary conditions at the four boundaries as a function of time starting with a uniform
mixture composition and a temperature spike at the center of the computational domain which serves
as an ignition source. Calculations have been performed with second and fourthorder accurate
discretizations of the spatial derivatives, and with exact and approximate factorization methods in
order to determine the influence of the spatial discretization and the factorization errors on the
numerical solution. Some sample results are shown in Figures 1 and 2 which show the temperature
(u) and fuel mass fraction (v) distributions at different times (t) and the errors that result from the
approximate factorization of the twodimensional operator into a sequence of two onedimensional
ones. Although the approximate factorization errors are secondorder in time, their largest magnitude
occurs at the steep fronts of both u and v.
Figure 1. u (top left), v (top right), and errors in u (bottom left) and v (bottom right) at t =10. (The
top figures were obtained with the implicit, linearized compact method presented here Dt=0.04,
Dx=Dy=0.5 and q=0.5. The bottom figures correspond to the difference between the results of a
secondorder accurate approximate factorization method in both space and time and the compact,
linearized technique presented here.)
Numerical experiments indicate that the exact factorization technique proposed here is able to accurately
resolve the steep gradients which occur in both the temperature and fuel mass fraction as the
combustion front propagates through the combustible mixture and approaches the boundaries of the
domain. The differences between the results obtained with second and fourthorder accurate
discretizations of the spatial derivatives are largest at the steep fronts; in fact, fourthorder accurate
discretizations provide steeper fronts than secondorder ones. Moreover, the iterative, exact
factorization method presented here converges in only two or three iterations and is very efficient
because the onedimensional finite difference operators provide block tridiagonal matrices and their
solution was obtained by means of LU decomposition.
Numerical experiments have also been performed in order to determine the effects of q, i.e., the
implicitness parameter, and the allocation of the reaction terms to the onedimensional operators
which result from the factorization. These experiments indicate that the factorization errors of the
fullyimplicit method, i.e., q=1, are secondorder in time which is the same order as that of the
CrankNicolson method, i.e., q=1/2; however, the discretization errors of the former, i.e., first
order accuracy in time, are larger than those of the latter, i.e., secondorder, and the fullyimplicit
method predicts a slightly faster flame speed than the CrankNicoloson technique when the
combustion front is sufficiently far from the boundaries.
The allocation of the reaction terms to the onedimensional operators results in errors which are also
secondorder in time and the differences between different allocations and the one which splits the
reaction terms in equal amounts amongst the different spatial operators are largest at the steep
combustion front. However, equal allocation of the reaction terms to the onedimensional operators is
computationally more demanding that those which assign the source terms to only onedimensional
operator because the Jacobian matrix of the reaction terms has to be evaluated for the two one
dimensional operators. Furthermore, for symmetric problems such as the one considered here, the
unequal allocation of the reaction terms to the onedimensional operators may result in solutions
which are not symmetric because the nonsymmetric distribution of the reaction terms increases the
splitting errors associated with the factorization of the twodimensional reactiondiffusion operator
into onedimensional ones, although, if the source terms are allocated to only onedimensional
operator, the other one is of the diffusion type.
Figure 2. u (top left), v (top right), and errors in u (bottom left) and v (bottom right) at t = 20. (The
top figures were obtained with the implicit, linearized compact method presented here Dt=0.04,
Dx=Dy=0.5 and q=0.5. The bottom figures correspond to the difference between the results of a
secondorder accurate approximate factorization method in both space and time and the compact,
linearized technique presented here.)
ACKNOWLEDGEMENTS
This work has been supported by Project No. PB941494 from the DGICYT of Spain.
NUMERICAL SIMULATION OF METHANE/AIR TURBULENT POOL FIRE
D. Morvan, B. Porterie, M. Larini and J.C. Loraud
IRPHE UMR CNRS 138 60 rue J. Curie Technopôle de Château Gombert
13453 Marseille cedex 13 France email: morvan@unimeca.univmrs.fr
IUSTI UMR CNRS 139 5 rue E. Fermi Technopôle de Château Gombert
13453 Marseille cedex 13 France email: berni@iusti.univmrs.fr
This paper reports numerical simulations of unconfined pool fires represented as a Methane/Air turbulent diffusion flame which develops from a porous burner. The modelling of turbulent combustion is solved using an eddydissipation combustion model coupled with a RNG ke statistical model for the turbulent flow. Radiation heat transfer and soot formation have been taken into account using P1approximation and transport submodels which permit to simulate main phenomena (nucleation, coagulation, surface growth ...) present in such system. The set of coupled transport equations is solved numericaly using a finite volume method, the velocitypressure coupling is treated with a projection technique.
INTRODUCTION
The progress performed for the modelling of turbulent combustion permit now to simulate with a more realistic manner pool fires which develop in industrial or natural environment. Field models represent a promissing way to improve knowledge of fundamental physical mechanisms which are associated with fire spread. They must permit for example to understand the interaction between gas flow dynamics (turbulence, hydrodynamics instabilities ...), heat flow (by convection, diffusion and radiation) and chemical kinetics (reaction rate, chemical species mixing ...) [1]. Because of high fluctuation levels observed in turbulent flame, the classical moments method generally used in turbulent flows modelling cannot be used for the calculation of chemical species in turbulent combustion. Various models based generally on physical considerations permit to represent multiple interactions between turbulent structures and reaction rate. The reaction rate for example could be evaluated from characteristic turbulent time scale (calculated with turbulent kinetic energy and its dissipation rate) and the fluctuation of fuel mass fraction mixture (Eddy Break Up model for premixed flames) or the average mass fraction of the more deficient chemical species (Eddy Dissipation Combustion model for diffusion flames). An other approach based on probability density function (Pdf) permits also to take into account of finite rate chemical kinetics, the average mass fractions of chemical species are calculated from simple quadrature between instantaneous species mass fractions exprimed as a function of instantaneous mixing fraction and the Pdf generally evaluated from average mixing fraction and its fluctuation [2,3].
Pool fires and jet flames constitute two regimes which limit the behaviour of diffusion flames. They can be distinguished from dimensionless heat flux Q' and the dependence of flame height:
where r_{Ą}, C_{p}, T_{Ą} design respectively the density, the specific heat and the temperature of the ambient air, and , g, d the energy release by the chemical reaction, the gravity magnitude and the burner diameter.
If the dimensionless flame height L/d is ploted as a function of dimensionless heat flux Q', we can distinguish three zones:
for Q’ <1 the flame height L varies inversely proportionally to d 2/3(pool fire regime)
for 1 < Q’ < 100 the flame height L is independant of burner diameter d
for 100 < Q’ the flame height varies linearly to the burner diameter d (jet flame regime)
Therefore pool fires regime can be characterized by a dependence of flame height with burner diameter and small dimensionless heat flux Q'. Because of large base dimensions and small fuel mass rate, wildeland fires and industrial fires are allways classified as pool fires.
MATHEMATICAL FORMULATION
Field modelling of pool fires is based on the resolution of conservation equations (mass, momentum, energy ...) inside the region including the fire and the solid oxidizer situated near this region. The buoyancy flow generated by hot gas expension due to heat release, is sufficiently important that we can suppose that the regime is locally turbulent. Average procedure, mass, momentum, species (fuel, oxidant, combustion products and inert gas) and energy conservation could be writen using classical Favre procedure [1]. The radiation heat transfer is treated using P1approximation, which consists to consider the gas mixture as a gray medium [4] with a mean absorption coefficient as a function of mole fraction of combustion products, soot volume fraction and mean temperature [5]. The turbulent flow is evaluated using RNG ke statistical model which permits to improve significantly the results obtained from classical high or low Reynolds number k e model for the description of turbulent flow including weak turbulent and recirculating flow regions [6]. Modelling of turbulent combustion consists to define a relation between reaction rate, turbulent variables (k, e...) and average or mean square fluctuations of mass fractions (chemical species, mixing fraction ...). The present calculation is based on Eddy Dissipation Combustion Model (EDCM) which represents an adaptation of the Eddy Break Up model for diffusion flames [7]. In this case the mean reaction rate is evaluated from a characteristic turbulent time scale evaluated from the turbulent kinetics energy and its dissipation rate, and the average mass fraction of the deficient chemical species (fuel, oxidizer or combustion products).
The present study includes the numerical simulation of an unconfined turbulent pool fire in a configuration represented in Figure 1, the dimension of the computational domain is 4 x 4 m. From numerical results the effects of input fuel mass rate upon the gaz flow, the dimension of the flame are studied.
ACKNOWLEDGMENT
The European Economic Commission is gratefully acknowledged for providing partial funding for this research in the frame of the EFAISTOS project.
REFERENCE
 Cox, G., Combustion fundamentals of fire. Academic Press, 1995.
 Borghi, R., Turbulent combustion modelling, Prog. Energy Combust. Sci., Vol. 14, pp 245292, 1988.
 Libby, P.A., Williams, F.A., Turbulent Reacting Flows, Academic Press,1993.
 Fusegi, T., Farouk, B., Laminar and turbulent natural convectionradiation interaction in a square enclosure filled with a nongray gas, Numerical Heat Transfer, Part A, Vol. 15, pp 303322, 1989.
 Kaplan, C.R., Shaddix, C.R., Smyth, K.C., Computations of enhanced soot production in timevarying CH4/Air diffusion flames, Combust. Flame, Vol. 106, pp 392405, 1996.
 Zijlema, M., Segal, A., Wesseling, P., Finite volume computation of incompressible turbulent flows in general coordinates on staggered grids, Int. J. Numerical Methods Fluids, Vol. 20, pp 621640.
 Magnussen, B.F., Hjertager, H., On mathematical modelling of turbulent combustion with special emphasis on soot formation and combustion, Proceedings 16^{th} Symposium (International) on Combustion, The Combustion Institute, pp 719729, 1976.
