Application of a flamelet-based CFD combustion model to the LES simulation of a diesel-like reacting spray

Application of a flamelet-based CFD combustion model to the LES simulation of a diesel-like reacting spray

Journal Pre-proof

Application of a flamelet-based CFD combustion model to the LES simulation of a diesel-like reacting spray ´ ´ J.M. Desantes, J.M. Garc´ıa-Oliver, R. Novella, E.J. Perez-S anchez PII: DOI: Reference:

S0045-7930(19)30377-9 https://doi.org/10.1016/j.compfluid.2019.104419 CAF 104419

To appear in:

Computers and Fluids

Received date: Revised date: Accepted date:

26 March 2019 5 November 2019 24 December 2019

´ ´ Please cite this article as: J.M. Desantes, J.M. Garc´ıa-Oliver, R. Novella, E.J. Perez-S anchez, Application of a flamelet-based CFD combustion model to the LES simulation of a diesel-like reacting spray, Computers and Fluids (2020), doi: https://doi.org/10.1016/j.compfluid.2019.104419

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2019 Published by Elsevier Ltd.

Highlights • Spray A from ECN is analysed with LES simulations and a flamelet model. • Results in very good agreement with experimental measurements. • Low fluctuations for the lift-off length due to the intense chemistry. • Ignition kernels appear spontaneously and detached from the main flame.

1

Application of a flamelet-based CFD combustion model to the LES simulation of a diesel-like reacting spray J.M. Desantes a , J.M. Garc´ıa-Oliver a , R. Novella a , E.J. P´erez-S´ anchez a

b

b,∗

CMT - Motores T´ ermicos, Universitat Polit` ecnica de Val` encia Edificio 6D, Camino de Vera s/n, 46022, Valencia, Spain Tel. (0034) 96 387 76 50 / Fax (0034) 96 387 76 59

Barcelona Supercomputing Center - Centro Nacional de Supercomputaci´ on Nexus II Building, C/ Jordi Girona 29, 08034, Barcelona, Spain Tel. (0034) 93 413 77 16 / Fax (0034) 93 413 77 21

Abstract Spray A from ECN, representative of diesel-like sprays, is modelled in the frame of Large-Eddy Simulations (LES) with a Dynamic Structure (DS) turbulence model in conjunction with an Unsteady Flamelet Progress Variable (UFPV) combustion model. In this work, the spray flow field is first calibrated under inert conditions against experimental data. In a second step, the reactive spray is simulated in order to describe the flame internal structure when varying ambient temperature. The model shows a good agreement with experimental results and describes the trends observed in flame global parameters, such as ignition delay (ID) and lift-off length (LOL). Low fluctuations are observed in LOL positioning revealing an intense chemical activity at the height of the base of the flame, which stabilizes the reaction in spite of the turbulent fluctuations. The analysis of the LES instantaneous fields shows how ignition kernels appear upstream of the base of the flame, clearly detached from the reaction zone, and grow and merge with the main flame in agreement with previous reported experimental and modelling results. The ambient temperature has a clear impact on the flame structure described by the model and the whole set of results reveal that in the frame of LES simulations the UFPV model is suitable for the calculation ∗ Corresponding

author. Email: [email protected]

Preprint submitted to Computers & Fluids

January 9, 2020

of diesel flames. Keywords: Large-Eddy Simulation, Spray A, non-premixed flames, chemical mechanism

1. Introduction The ever-increasing relevance of the transport and energy sectors in our society has lead to the optimization of combustion devices in order to increase their efficiency and reduce pollutant emissions. In particular, understanding the com5

plex phenomena involved in diesel-like reacting sprays and their interaction is a challenging field in the research of partially premixed and non-premixed combustion. These phenomena could be summarized in atomization and break-up, evaporation, mixing, chemical oxidation and spray-wall interaction occurring in a high Reynolds turbulent flow [1].

10

The complete resolution of the whole physical and chemical processes developing in so different spatial and temporal length scales leads to an unaffordable computational cost and, consequently, different hypothesis have to be introduced in order to derive simplified models. However, these assumptions entail new uncertainties that have to be verified. Unfortunately, measuring in an en-

15

gine is a difficult task that, in general speaking, only provides global or integral variables, insufficient to give an exact picture of the whole combustion process and validate the models. In this context, the Engine Combustion Network (ECN) [2] has suggested a set of representative experiments to be carried out in special combustion

20

chambers, constant-volume pre-burn (CVP) combustion vessels and constantpressure flow (CPF) rigs [3, 4], that discard many uncertainties and allow to measure with the most advanced experimental techniques. The empirical observations are complemented with Computational Fluid Dynamics (CFD) simulations. The participation of a wide sector of the researcher community allows

25

to diffuse and improve the state of the art. In line with this goal, this work deals with simulating the well-known spray

3

A from ECN, a single nozzle spray with boundary conditions corresponding to modern diesel-like sprays. Liquid length, vapour penetration and some relevant spatial fields together with ignition delay (ID), lift-off length (LOL) are mea30

sured in these combustion chambers. A CPF rig is available at CMT-Motores T´ermicos whose experimental results are used along this paper [5, 6] unless otherwise stated. Regarding diesel spray modelling an extensive literature may be found for inert and reactive conditions. As mentioned previously, the great variety of

35

physical and chemical phenomena occurring at different scales in the flow has given rise to the formulation of many different models aiming to describe these phenomena with different levels of accuracy and computational cost. For years, the Reynolds Averaged Navier-Stokes (RANS) equations have been widely used due to their relative reduced computational cost compared to

40

other approaches. Excellent results have been reported in the frame of combustion modelling by means of RANS simulations [7, 8, 9, 10]. Notwithstanding, and despite the positive capabilities of the RANS approach, as the whole range of scales are modelled their hypothesis may not be completely fulfilled. In order to solve more accurately the flow, Large-Eddy Simulations (LES)

45

have gained attention during the last years since, in spite of their higher computational cost, the large eddies are solved and only the smallest scales, that tend to show statistical isotropy and universality [11], have to be modelled leading to the development of more accurate models [12]. Notwithstanding, chemical reactions occur when the species are mixed at

50

molecular level implying that combustion takes place at or even below the smallest scales of the flow. Hence, the LES simulation does not solve the chemical source term, which is completely modelled [13, 14], and this leads to directly extend the RANS combustion models to the LES approach. However, it is expected to obtain more accurate results for LES since the shape of the filtered

55

probability density function (FPDF) has less impact [14] and the integral scales are better predicted. These advantages are conjugated with the capability to reproduce intermittency [15, 16]. LES calculations have been successfully applied 4

to combustion simulation [17, 18, 19, 20]. Regarding the LES turbulence models, some of them are based on the Boussi60

nesq hypothesis as the Smagorinsky model, which requires adjusting ad hoc the value of its constant CS . This shortcoming can be palliated with the dynamic model [21]. More refined approaches can be obtained when transporting the sub-grid kinetic energy ksgs as in the One Equation Eddy model (OEE) [22], which allows to use coarser meshes which are more compatible with Lagrangian

65

droplet models, conventionally used to model the liquid phase in spray simulations. In addition, the presence of two phases in diesel sprays adds new terms to the sub-grid scale modelling and complicates the turbulence models [23]. Finally, in other approaches like the Dynamic Structure model (DS) [24], which is not based on the Boussinesq hypothesis, the residual stress fluxes are

70

considered proportional to the sub-grid turbulent kinetic energy ksgs by means of a coefficient tensor and ksgs is transported. This model fulfils some important properties that are desirable for any turbulence model, namely, solvability, scaling, frame invariance and realizability [24]. Due to its capabilities to perform accurate simulations and its formulation for multiphase simulations [25, 26], it

75

has been adopted in the current work. Solving a diesel spray simulation requires a combustion model with the ability to manage complex chemical schemes while retaining turbulence-chemistry interaction (TCI). In some works the diesel spray has been modelled with Perfect Stirred Reactors (PSR) [27, 28] but, as Bhattacharjee et al. [29] point out,

80

the lack of TCI makes this model not suitable for this problem especially for conditions of low reactivity. It is worth mentioning the previous mentioned works were carried out in the frame of RANS simulations, however, it is expected that when calculating with LES simulations the impact of not considering TCI, even not being a correct modelling approach since chemistry occurs at the small-

85

est scales of the flow which are not solved in LES, may lead to more accurate results [17]. Other models like the Representative Interactive Flamelet (RIF) model can provide good results [30, 31, 32] but lead to a non-local description of the flame and the need of injecting a high number of flamelets, that increase 5

the computational cost in order to obtain accurate results. Using similar con90

cepts the Flamelet Generated Manifold (FGM) has been demonstrated to be a powerful model that captures many phenomena of the diesel flame [33, 20]. Finally, transporting probability density functions (PDFs) has been shown to be a fruitful approach to describe combustion, for which very accurate results for ID and LOL as well as flame structure have been provided [29, 9, 34]. How-

95

ever, its very high computational cost makes this model difficultly affordable for practical engine simulations. Regarding the physical phenomena observed in the diesel flame a great effort has been devoted to study combustion development at the LOL. This distance delimits the height at which important amounts of heat are released and provides

100

an indirect measure of the soot being formed, since it determines the maximum mixture fraction that reacts in the spray [35]. In this region the flame stabilizes and the mechanism of how this occurs is still a matter of intense study. Experimentally it has been observed how hot pockets of partially premixed fuel-air mixtures ignite spontaneously upstream and detached from the flame [36, 37]

105

and this observation has been confirmed by means of LES simulations [17, 20]. This is an important result since it is an indirect evidence that one of the flame stabilization mechanisms is auto-ignition although some authors specify that the predominance of auto-ignition or other stabilization mechanisms, such as flame propagation, depends on the boundary conditions too, being more important

110

flame propagation for high ambient temperatures [38]. Due to its interest this aspect will be covered in this paper too. Considering previous picture for the different models available in the literature and their capabilities, the flamelet-based model appears as one of the most promising alternatives [39, 40, 41, 42, 43, 8] with a remarkable balance between

115

accuracy and computational cost. This model describes the turbulent flame as an ensemble of strained laminar flames called flamelets [44]. The flamelet equations may be rewritten in the mixture fraction space yielding a system of one-dimensional transport equations [45, 46]. In LES framework, the solutions provided by this system are integrated by means of FPDFs in order to account 6

120

for the TCI. Depending on the approach, the flamelet equations may be solved during the CFD calculation, as in the Representative Interactive Flamelet model (RIF) [30], or, on the contrary, the flamelet solutions may be pre-tabulated, as in the Flamelet Generated Manifold (FGM) [39] or the Flame Prolongation of ILDM

125

(FPI) [47], which are based on the Intrinsic Low Dimensional Manifold (ILDM) [48]. The FPI approach is adopted in this work. However, despite the advantages provided by these models the computational cost may increase exponentially when dealing with diesel engine simulations where the boundary conditions may span over wide ranges, complex mechanisms

130

are required and a local description of the flame is necessary [49, 50]. In this context, a simplified approach called Approximated Diffusion Flamelets (ADF) was suggested in order to reduce drastically the computational cost while retaining the ability to manage complex chemistry [49]. This approach is adopted for the current work given the excellent reported results provided in the literature

135

[51, 43, 50, 10]. With this work many results found by means of experiments and modelling techniques [52, 33, 53, 36, 37] are intended to be corroborated with the application of the Dynamic Structure model used in conjunction with the ADF flamelet model which is used in the context of LES simulations. To demon-

140

strate the suitability of such combination of modelling approaches applied to diesel sprays simulations is another target which is intended to be covered with the perspective of giving validity to future numerical experiments where more complex phenomena will be incorporated. Hence, the objective of this work is the analysis of the diesel flame structure, with special attention to the behaviour

145

at the base of the flame where the LOL is positioned, and how such structure is influenced by the change of the boundary conditions, in particular the air temperature. For this purpose, spray A from ECN [2] has been simulated since there exists a vast set of experimental results that allows to contrast many of the conclusions extracted along the work.

150

The paper is structured with the description of the problem, followed by 7

a section where the numerical algorithms and mesh are explained. This information is complemented with the description of the DDM (Discrete Droplet Method) and combustion models. Subsequently, a first validation in inert conditions is carried out which is closed with an initial assessing of the reacting 155

flame. The capability of the model is evaluated with the analysis of the global parameters, ID and LOL, with special attention to the last one. The instantaneous fields are described by means of spatial and Z-T maps, first for the nominal case, and then for the rest of boundary conditions. This is followed by an analysis of the flame at the LOL in order to shed light in the flame stabiliza-

160

tion mechanism. Previous analysis is complemented with the description of the averaged fields before giving the conclusions of the work which close the paper.

2. Problem description The objective of this work is the simulation of spray A in the LES framework by means of the flamelet concept. The boundary conditions for spray A are re165

produced in table 1 corresponding to a temperature parametric sweep. Nominal injector diameter is 90 µm, with nozzle code 210675 [2] and discharge coefficient equal to 0.9 [54]. The fuel is n-dodecane which is used as a diesel surrogate. A long injection mass rate (4 ms) with an injection pressure pinj of 150 MPa is imposed [55] where the fuel temperature is assumed equal to 363 K. In table

170

1 the values for the stoichiometric and saturation mixture fractions, Zst and Zs , respectively, are included. The value for Zs corresponds to the maximum mixture fraction value for which fuel does not condensate in the air-fuel mixture and is calculated from a mixing-controlled approach [56]. Table 1: Definition of spray A parametric study. XO2

Tamb (K)

ρamb (kg/m3 )

pamb (MPa)

pinj (MPa)

Zst

Zs

0.15 0.15 0.15 0.15

750 800 850 900

22.8 22.8 22.8 22.8

4.97 5.3 5.63 5.96

150 150 150 150

0.046 0.046 0.046 0.046

0.251 0.278 0.303 0.326

The CFD simulation is carried out in the open CFD platform OpenFOAM 8

175

environment [57] with an in-house developed code. The mesh is a cylinder composed of approximately 3.6 Million cells with a definition in terms of cell size and distribution based upon the spray geometry. The cylinder has a radius of 23.5 mm and height of 108 mm with open boundary conditions for all the faces except the base of the cylinder, where the injector is positioned, which is

180

a wall. The mesh has an O-grid structure and is composed of an internal prism with its axis coincident with the cylinder one. The length of the side of the base of such prism measures 1.688 mm and is composed of square cells with a constant size of 62.5 µm, following the recommendation given in the literature for the use of Lagrangian droplet models [58, 59]. The rest of the mesh follows

185

a discretization in cylindrical coordinates with 108, 108 and 292 cells in radial, azimuthal and axial directions, respectively, and an expansion ratio between cells of 1.015 for the radius and 1.01 for the axial distance. Different views and cuts of the mesh are shown in figures 1 and 2.

Figure 1: 3D view of the mesh and the domain.

Regarding the droplets modelling, the carrier phase is fully coupled to the 190

Lagrangian phase applying a two-way coupling. In addition, a source term for the turbulent kinetic energy is contemplated too [25, 60]. A DS model has been adopted to solve the residual stresses in the diesel spray simulation in the formulation described in [25, 60] and the implementation given in [26], which accounts, as said before, for the sub-grid interactions between the gas phase and 9

1,688

47

108

Figure 2: Different cuts of the mesh used for the calculations. Top figure shows a meridian plane. Bottom figure shows a perpendicular cut to the cylinder axis with a zoom at the inner prism. Dimensions in mm.

195

the fuel droplets. For these droplets a Lagrangian approach has been used with the Reitz model [61] that models the Kelvin-Helmholtz and Rayleigh-Taylor instabilities. The constants applied to this model for the current simulations are B0 = 0.61, B1 = 10, Cτ = 1 and CRT = 0.1 with the Weber number limit fixed at 6.

200

Moreover, the half angle of the spray cone is set equal to 12o . The droplets are injected following a Rosin-Rammler distribution where the droplet diameter ranges from 10 µm to the diameter of the nozzle and the distribution is defined with an exponent n = 3. In this way, a population of droplets is injected into the domain and no atomization model is used since there is no liquid vein to

205

break. The evaporation is solved with a Ranz-Marshall correlation but since there exists evaporation, the Nusselt number is recomputed by means of Bird’s correction. Finally, collision is treated with the stochastic O’Rourke model [62]. The combustion model is based on the flamelet concept where the set of

10

laminar flamelets are solved in the mixture fraction Z space. According to 210

the theory devised by Peters [45, 46], for each flamelet a system of partial differential equations (PDEs), that accounts for the transport phenomena, is solved in conjunction with the ordinary differential equations (ODE) system for the chemical source terms. The flamelet equation reads for each species

∂Yk χ ∂ 2 Yk + ω˙ k = ∂t 2 ∂Z 2

k = 1, . . . , N

(1)

where Yk is the species mass fraction for species k and N is the total number 215

of species. The scalar dissipation rate χ measures the diffusivity in the Z space and is prescribed following the profile

χ(a, Z) =

F (Z) a 2 Z exp[−2(erfc−1 (2Z/Zs ))2 ] = χst π s F (Zst )

(2)

where subscript st refers to stoichiometric values and a is the strain rate [46]. The set of equations given by the system (1) together with the chemical ODE system may involve a very high computational cost when solving complex 220

chemical mechanisms and a wide set of boundary conditions, as it is the case of diesel engines simulations. Based on these limitations the ADF model was suggested in order to deal with complex mechanisms in reduced computational times when solving the flamelet equation while still retaining the flame structure [49]. In this approach

225

only the equation for the progress variable Yc is solved, where Yc is defined as a linear combination of species mass fractions. The chemical source terms as well as the relationships between the species and Yc come from solving a set of homogeneous reactors (HRs). This leads to solve first the corresponding set of HRs, which are not computationally expensive, where the initial conditions

230

come from the adiabatic mixing between air and fuel with the same boundary conditions than the flamelet [56]. Such HRs are solved imposing adiabatic

11

boundary conditions and constant pressure [49]. In this work this has been done with Chemkin software [63]. Then the flamelet equation to be integrated reads

∂Yc χ(a, Z) ∂ 2 Yc + ω˙ cHR (Z, Yc ) = ∂t 2 ∂Z 2

(3)

where the dependencies are explicitly written. Equation 3 is solved for dif235

ferent values of strain rate. This approach has been adopted to solve spray A in this work and its implementation has been carried out following [64]. The steady solutions are first computed with a Newton-Rapson algorithm and then the unsteady solutions are calculated with an implicit method and an adaptive time step algorithm based on the chemical source term. Second order

240

is used for mixture fraction derivatives while first order for temporal derivatives. Once the flamelet equations are computed, a database in the form ψ = ψ(Z, χst , t) is available where ψ is any reactive scalar. TCI, necessary to model the sub-grid effects of turbulence, is accounted for by means of FPDFs. Follow-

245

ing a similar approach to that used in other works [65, 51, 43], a beta function g 002 ) where Z 002 e Zg is used in the mixture fraction direction [46], defined by (Z, sgs

sgs

corresponds to the sub-grid mixture fraction fluctuations, while for the rest of variables δ-functions are assumed. Additionally, statistical independence is hypothesized leading to

002 f , e e Z, e Zg ψ( st t ) = sgs , χ

Z∞ ZZsZ∞ 0

0

0

002 e Zg ψ(Z, χst , t) δ(t − e t) PZ (Z, Z, f st ) dt dZ dχst (4) sgs ) δ(χst − χ

which becomes

002 f , e e Z, e Zg ψ( st t ) = sgs , χ

250

Z

0

Zs

002 e Zg ψ(Z, χf t ) PZ (Z, Z, st , e sgs ) dZ

(5)

002 f , Y e Z, e Zg ec ) providing a turbuEquation (5) is reparametrized as ψe = ψ( st sgs , χ

12

lent flame manifold. In the same way, χ profile is integrated yielding a relationship in the form

002 e g χ e = χf st J(Z, Zsgs )

(6)

002 e 17 for Zg In this work, around 32 values have been tabulated for Z, sgs , 35

255

e for χf st , depending on the boundary conditions, and 51 for Yc with a parabolic distribution in order to retain accurately auto-ignition.

Continuity, Navier-Stokes and energy equations are solved together with the

sub-grid kinetic energy. Additionally, the most important species (H, OH, CO, CO2 , H2 O, C12 H26 , CH2 O, C2 H2 , C7 H14 , H2 , O2 , N2 ) are transported [66] and

260

f the set of species chemical source terms ω ˙ k is computed from the turbulent 002 have to e and Zg flame manifold. For such purpose, transport equations for Z sgs

be solved which read

e e ∂(ρ Z) ∂(ρ uei Z) ∂ + = ∂t ∂xi ∂xi 002 002 ∂(ρ uei Zg ∂(ρ Zg ∂ sgs ) sgs ) + = ∂t ∂xi ∂xi

∂ Ze ρ Def f ∂xi

002 ∂ Zg sgs ρ Def f ∂xi

!

!

+ SZ

+ 2ρ Dsgs

(7)

e ∂ Ze ∂ Z − ρ χg sgs (8) ∂xi ∂xi

Def f is the effective mass diffusion coefficient which is computed as the sum of the laminar and sub-grid (Dsgs ) contributions. SZ is the source term for

265

mixture fraction caused by evaporation. 002 equation the sub-grid dissipation term χ In the Zg g sgs is modelled by [51, 43] sgs

χg sgs = Cχ Dsgs

002 Zg sgs ∆2

(9)

As said before, Dsgs is the sub-grid mass diffusivity defined from νsgs with a Schmidt number 0.7 [25], Cχ is a constant to be adjusted and ∆ is the size of the 13

e2 cell. The total scalar dissipation rate χ e = χg sgs + 2D|∇Z| , with D the laminar 002 ) by means of equation 6. e g mass diffusivity, is related to χf st by the pair (Z, Z sgs

270

From the transported species mass fractions Yec is obtained, which is defined 002 f , Y e g ec ) the ] as Yec = Yg CO + YCO2 [49]. With the input parameters (Z, Zsgs , χ st

source term for Yec , ∂ Yec /∂t, obtained from solving the flamelets and integrating,

is retrieved and used to obtain the Yec at next time step due to chemical effects according to the following expression

275

∂ Yec e g 002 , χ e Yec (t + δτ ) = Yec (t) + (Z, Zsgs f st , Yc (t)) δτ ∂t

(10)

where the different values are taken in each cell of the domain. From 002 f , Y e Zg ec (t + δτ )) and the corresponding tabulation, the species mass (Z, st sgs , χ

fk fractions Y

tab

are also retrieved from the look-up table, enabling the calcula-

f tion of the species source terms ω ˙ k according to the following equation. f ω ˙k =

280

fk Y

tab

002 f , Y e Zg ec (t + δτ )) − Y fk (Z, st sgs , χ δτ

(11)

fk is the mass fraction for species k at the cell. δτ is the time step where Y

for advancing in the chemical manifold which is taken equal to the CFD time step, fixed to 0.02 µs. This approach, even being more expensive, has the advantage that it eases handling additional effects like evaporation or pollutant emissions compared to other models where only the mixing and progress variables are transported and

285

the species are retrieved from the manifold [66].

3. Numerical modelling The flow is solved by means of the finite volume method. To solve the transport equations the PISO (Pressure Implicit with Splitting of Operators)

14

algorithm is applied where the temporal derivatives, the Laplacian and diver290

gence terms are evaluated with second order schemes. For the current calculations, the dodecane oxidation is described by the chemical scheme developed by Narayanaswamy et al. [67] which comprises 255 species and 2289 reactions, which is widely extended in the literature [68, 69].

4. Results and discussion 295

The analysis is divided in a first section where the model is validated in terms of the nominal inert and reactive spray and a second part devoted to the description of the reactive spray and how is influenced by the boundary conditions. 4.1. Validation of the model

300

The spray is first calibrated for the inert nominal condition by means of modelling and experimental data comparison. As a preliminary step, the vapour penetration and the liquid length for the simulated and the experimental results [6] are shown in figure 3. The liquid length is defined as the distance to the nozzle where 95% of the injected liquid is found and the vapour penetration is

305

given by the maximum distance from the nozzle outlet to where mixture fraction is 0.001 [2]. Although the liquid length is slightly overestimated an excellent agreement is observed for the spray vapour penetration since simulated and experimental curves fall very close.

310

The validation of the simulation is followedq by the comparison of the mean g 002 i or Z e its standard deviation, hZ mixture fraction, hZi, rms , and the mean axial velocity he ui normalized by the velocity at the exit of the nozzle. The symbol h i denotes averaged value. These profiles obtained from the simulation are averaged in time and azimuthal directions and are compared with experimen-

315

tal results in figure 4. The measurements were obtained by means of Rayleigh imaging for mixture fraction and Particle Image Velocimetry (PIV) for velocity

15

100

experiment simulation

S [mm]

80 60 40 20 0

0

1

2 t [ms]

3

4

Figure 3: Vapour penetration and liquid length for experimental (black lines) and simulated (blue lines) inert nominal condition. For experiments, uncertainty of measurements is delimited with shadows.

[70, 71, 72, 53]. PIV measurements were carried out in the interval that extends from 40 to 100deq , approximately. As the measurements and the simulations

320

have been carried out with different nozzles the distances are normalized by the p equivalent diameter deq = d0 ρf /ρa where d0 is the nozzle diameter and ρf , ρa are the fuel and air densities, respectively. The radial cuts at 50 and 90 deq correspond to axial distances of 25 and 45 mm, approximately. e and he An excellent agreement is observed for hZi ui. For Zrms , an analysis

of the constant Cχ shows that it has no strong impact on the results since such 325

constant only affects the sub-grid component. Notwithstanding, it is observed a slight decrease of Zrms when increasing the Cχ constant and, finally, a value of Cχ = 5 has been chosen to calibrate the model. With this value the experimental profiles for Zrms are reasonably well-captured on the axis and radial cuts. Moreover, the energetic spectral content of turbulence is gathered in figure

330

5 for two points on the axis and other radially displaced for the inert nominal case. The energy spectrum is obtained in the frequency domain, following [73], from the axial velocity signal registered at a given point during a temporal window where the turbulent statistics stabilize. Then the velocity autocorrelation function is computed and the Fourier transform is applied for this function

16

0.2

0.6

Zcl inert exp Ucl inert exp Zcl inert sim Ucl inert sim

0.15

experiment simulated

0.15 x/deq = 50

0.45

0.3

0.05

0.15

Z [−]

0.1

Ucl/U0 [−]

Zcl [−]

0.1 x/deq = 90 0.05

0

40

60

80 x/deq [−]

100

0 −15

0

−10

−5

0 r/deq [−]

5

10

0.05

0.07 0.06

experiment Cχ = 3

0.05

Cχ = 5 Cχ = 7

0.04

0.04 Zrms [−]

Zrms [−]

120

0.03

15

experiment Cχ = 3

x/deq = 50

Cχ = 5 Cχ = 7

0.03 0.02

0.02

0.01

0.01

x/d

eq

0 0

20

40

60 80 x/deq [−]

100

120

0 −15

−10

−5

0 r/deq [−]

5

= 90 10

Figure 4: Spray validation for very advanced time steps for inert nominal condition. Top e and normalized he e radial profiles at 50 and 90 deq . left: hZi ui along the axis. Top right: hZi Bottom left: Zrms along the axis for different Cχ values. Bottom right: Zrms radial profiles at 50 and 90 deq for different Cχ values. Experimental uncertainties delimited with shadowed regions.

335

providing the energy spectrum. From figure 5 is clearly seen that there exists a range of frequencies where the spectrum is almost flat corresponding to the integral scales of the flow. In spite of the noise in the signal, the energy spectrum is observed to fall in the high frequency range with a slope close to the value -5/3 predicted by the classic

340

turbulence theory, which corresponds to the inertial sub-range [11, 12]. This range is positioned in frequencies between 104 and 105 Hz. Besides, figure 5 shows that when moving away from the nozzle on the axis or radially away from the axis, the inertial sub-range tends to shift towards lower frequencies, which implies that the flow dynamics become slower as a result of the reduction of the 17

15

2

10

2 2

log Ey [m /s ]

0

10

−2

10

−4

10

2

10

3

10

4

10 log f [Hz]

5

10

6

10

2

10

2 2

log Ey [m /s ]

0

10

−2

10

−4

10

2

10

3

10

4

10 log f [Hz]

5

10

6

10

2

10

2 2

log Ey [m /s ]

0

10

−2

10

−4

10

2

10

3

10

4

10 log f [Hz]

5

10

6

10

Figure 5: Energy spectrum in the frequency domain for points on the axis at a distance of 20 (top) and 30 mm (center) and a 2 mm radially displaced point positioned 30 mm from the nozzle (bottom). For reference a straight line with slope -5/3 is included too. Results for inert nominal case.

18

345

velocities. In general, this range intends to be partially or completely solved with LES simulations. It is worth mentioning that it is probable that the current simulations do not solve completely the inertial sub-range because of the large cell sizes that the DDM formulation forces in order to fulfil its hypotheses. As pre-

350

viously said, in this work a minimum cell size of 62.5 µm has been imposed following the results obtained in the literature for diesel sprays when using the DDM approach [58, 59]. Finally, to complement previous validation of the mesh and configuration of the model, figure 6 shows the fraction of the sub-grid turbulent kinetic energy

355

in the resolved motions for the inert nominal case at an advanced instant, that is,

M=

hksgs i hksgs i + hksolved i

(12)

where ksgs is the sub-grid turbulent kinetic energy, which is obtained solving its transport equation, and ksolved is the solved turbulent kinetic energy computed as

ksolved = 360

1 [(u1 − hu1 i)2 + (u2 − hu2 i)2 + (u3 − hu3 i)2 ] 2

(13)

where in all previous expressions h i denotes an average performed in azimuthal and temporal directions and ui with i = 1, 2, 3 are the components of velocity. According to Pope’s criterion [13], it is desirable that M is lower than 0.2 at all the points of the mesh. In figure 6, it is observed that, in the vicinity of the stoichiometric mixture fraction, M values around 0.2 and 0.3, which even

365

being a little higher than Pope’s criterion, is still reasonable and acceptable. It is important to note that the Lagrangian model adopted for the droplets prevents from reducing the cell size at will since then the gas volume fraction at each cell would be too low and the hypotheses of the DDM would not be fulfilled [74]. Far away of the stoichiometric mixture fraction at the head of the

370

spray the values reduce to 0.1 and very far from the spray, where velocities are 19

r [m]

0.02

0.4

0.01

0.2

0

0

0.02

0.04 0.06 x [m]

0.08

0.1

0

Figure 6: Map for M (see equation 12) for the inert nominal case at an advanced instant. Black line represents the stoichiometric mixture fraction level curve.

0.25 experiment simulation

experiment simulation

80

0.2

60

0.15

Ucl/U0 [−]

S [mm]

100

40 20

0.1

0.05

0 0

1

2 t [ms]

3

4

0

40

60

80 x/deq [−]

100

120

Figure 7: Left figure: tip penetration and liquid length for experiment and simulation. For experiments, uncertainty of measurements is delimited with shadows. Right figure: normalized he ui on the centreline. Results for reactive nominal condition.

extremely small, M increases again. In this region the mesh is coarser in order to reduce the computational cost what explains the increase of M , nevertheless, it is deemed that such increase has not an important impact on the results due to the low velocities found in this region. 375

The calibration of the inert spray is considered quite satisfactory and it encourages the validation of the spray in reactive conditions. Again this is done by the comparison of modelling and experimental data for the reactive nominal case. Figure 7 gathers the comparison of the tip penetration and the liquid length as well as the axial velocity on the axis normalized by the velocity at

380

the exit of the nozzle. The tip penetration for the reactive case is defined in a similar way that the vapour penetration, that is, the maximum distance from the nozzle outlet to where mixture fraction is 0.001. Figure 7 shows that there exists an excellent agreement for the tip pene-

20

tration. Regarding the liquid length is slightly overestimated as in the inert 385

condition. However, as the liquid region is spatially isolated from the flame due to spray A boundary conditions it is expected that this overestimation does not influence subsequent results. Finally, the axial velocity on the axis from the simulation shows an acceptable correspondence with experimental results. The validation of the model is closed with a qualitative comparison of formalde-

390

hyde (CH2 O) and hydroxide (OH) fields with measured Laser Induced Fluorescence (LIF) data from [53] at advanced instants for the reactive nominal case. Experimental results are averaged in time and, consequently, are compared with averaged LES simulations results. Additionally, for the sake of completeness, instantaneous snapshots for the LES simulation are included too. Since the nozzle

395

diameter in experiments is slightly different from that used for the calculation, the distances are again normalized by the equivalent diameter. These results are gathered in figure 8.

21

25 r/deq [−]

r/deq [−]

25

0

−25

−25 20

40

60 x/deq [−]

80

100

0

r/deq [−]

r/deq [−]

0

25

0

0

20

40

60 x/deq [−]

80

20

40

60 x/deq [−]

80

100

60 x/deq [−]

80

100

80

100

25

0

100

0

20

40

25 r/deq [−]

25 r/deq [−]

0

0

−25

0

−25 0

20

40

60 x/deq [−]

80

100

0

20

40

60 x/deq [−]

Figure 8: Comparison of the CH2 O and OH fields for the nominal case. First row corresponds to experimental results [53], second row is for LES averaged solutions and third row is for LES snapshots. Left column shows 355 nm LIF signal for the experiment and the CH2 O field for the simulation while right column is derived from OH LIF signal for the experiment and shows the OH field for the simulation. LOL values for experiment (white dashed) and simulation (red dashed) are shown and the stoichiometric level curve (solid white line) is included for the modelled results.

CH2 O is detected in the vicinity of LOL in both experimental and simulated results. In the first case, the signal has been saturated downstream of 400

50deq in order to better visualise the CH2 O field since, as discussed in [53], this downstream signal is probably due to the interference of polycyclic aromatic hydrocarbons (PAHs). LES simulation predicts a narrow region of CH2 O production in the zone of rich mixtures close to the LOL and positioned in the axial distance of 30-50 deq in agreement with experimental results.

405

It is observed that in the experiment the CH2 O field seems to extend upstream of the measured LOL while in the simulation it appears downstream of 22

the LOL obtained from CFD. It is difficult, however, to extract any conclusion from this fact since LOL definitions for experiments and simulations are not exactly the same so this comparison is intended to be only qualitative and, hence, 410

some discrepancies are expected to arise. Regarding the OH field the simulation predicts two regions of OH production in a similar way than the experiment although in the experiment OH spreads over a wider region. In both cases OH is produced downstream of the corresponding LOL position. In the experiment the laser sheet extends until 92deq ,

415

approximately, and no further information is available. The results obtained so far show that there exists a good agreement between the LES simulations in inert and reactive conditions and the experiments and, hence, this enables to carry out the analysis of the reactive spray. 4.2. Analysis of the reactive flame

420

The analysis of the reactive spray is carried out in terms of flame metrics, namely ID and LOL, and reactive scalar fields in spatial and mixture fractiontemperature (Z-T) space representations for the different boundary conditions. 4.2.1. Global parameters First, the ID and LOL for the temperature sweep are shown in figure 9. The

425

ID is defined as the time spent from start of injection (SOI) until the maximum Favre filtered ambient temperature increases 400 K [1]. For the LOL a first definition, following the ECN criterion, has been adopted [2], which defines this value as the minimum axial distance from the nozzle to the surface level given by the averaged Favre filtered OH mass fraction corresponding to the 14 % of

430

the maximum value reached in this field. This LOL value is labelled as ‘LES OH aver.’ in figure 9. However, a second definition for the LOL has been considered which is computed as the average of the minimum axial distance from the nozzle to the surface level of the Favre instantaneous OH mass fraction corresponding to the 14 % of the maximum value reached in this field. This value is labelled

23

435

as ‘LES OH instant.’ in figure 9. This will make possible to evaluate the impact of the order of the operations in defining the LOL. Representative values for Yg OH are found in the vicinity of the stoichiometric

mixture fractions where the temperature is maximum. However, during the

440

injection rate, which lasts for 4 ms, the surface levels for Ze = Zst do not stabilize completely making difficult to obtain representative averaged Yg OH fields and,

particularly, the maximum value. Consequently, as the first methodology (ECN criterion) first averages the Yg OH field this may cause some uncertainties in the

LOL value and, in consequence, to avoid these shortcomings a second definition for the LOL has been suggested.

3

experiment model

experiment LES OH aver. LES OH instant.

50

2.5

40 LOL [mm]

ID [ms]

2 1.5 1

30 20

0.5

10 0 750

800

850

750

900

800

850

900

T [K]

T [K]

Figure 9: Ignition delay (left) and lift-off length (LOL) (right) for the temperature parametric sweep. Error bars delimit the standard deviation of the variable.

445

As observed in figure 9 there exists an excellent agreement between ID results for all the temperatures except for the lowest one, for which the LES simulation predicts a too fast ignition. Regarding the LOL, both criteria provide similar results although a slight deviation is accentuated when decreasing the reactivity. Moreover, both definitions give short LOL values although they reproduce the

450

expected trends. The LOL values for the second definition (‘LES OH instant.’) seem to show a clearer parallelism with regard to the experimental curve than the first one (‘LES OH aver.’). In addition, the standard deviation has been computed for the second defini-

24

tion and included in the figure by means of error bars. This fluctuation is com455

puted from the resolved motions and does not consider a sub-grid component. However, it is thought that this fact does not modify subsequent conclusions. It is observed that the fluctuation level is not very high in agreement with the experimental results which are of the same order of magnitude (±1 mm) than those corresponding to the simulations.

460

In order to visualize how the LOL evolves in time, figure 10 shows the temporal evolution of the LOL computed with the instantaneous Yg OH field. The

LOL value computed with the first criterion is included too. 60

750 K 800 K 850 K 900 K

50

LOL [mm]

40 30 20 10 0

0

1

2 t [ms]

3

4

Figure 10: Instantaneous LOL evolutions for the different ambient temperature cases (solid line) and LOL defined with ECN criterion (dashed line).

From figure 10 it arises that the first ignition kernels appear at positions downstream the distance where LOL finally stabilizes since LOL position moves 465

upstream during the first steps of ignition, as reported in other works [1, 29]. Once the LOL stabilizes its fluctuations are low (±1 mm) and do not depend clearly on the ambient temperature. Moreover, LOL position signal shows rapid and strong fluctuations, similar to discontinuities, when it stabilizes, e.g. at 2.7 ms for 900 K or 3.5 ms for

470

800 K. This is not due to a rapid recession of the base of the flame but to the apparition of ignition kernels that spontaneously start burning upstream and detached from the base flame as will be explained in detail in section 4.2.2.

25

Therefore, the base of the flame, defined by the LOL position, is anchored in the vicinity of a fixed point and does not experiment high fluctuations implying 475

that the diesel spray generates a stable and vigorous flame. An intense partially premixed flame is established in the vicinity of the LOL where reactants rapidly burn into intermediate products releasing important amounts of heat [52]. This entails that the chemical activity is weak upstream of the LOL but it is triggered at the LOL distance. The intensity of chemistry is so high that it is not

480

strongly affected by the local flow conditions and then LOL shows low fluctuations as observed in figures 9 and 10. A deeper insight into the occurrence of the mentioned ignition kernels at the LOL distance is developed in next section for the different boundary conditions. 4.2.2. Instantaneous fields

485

This section shows the instantaneous fields for some important reactive scalars in both spatial and Z-T representations for advanced instants and an analysis of the formation of the ignition kernels at the base of the flame. g Figure 11 gathers the temperature, Y^ CH2 O and YOH , which are tracers of the

low and high temperature combustion respectively, for instants 3500, 3800 and 490

4000 µs corresponding to the nominal condition.

g Figure 11: Instantaneous fields for Te (first row), Y^ CH2 O (second row) and YOH (third row) e = Zst is for the instants 3500 (left), 3800 (center) and 4000 µs (right). The level curve Z included in white line. LOL values for experiment (white dashed) and simulation (red dashed) are shown too. Results for nominal case.

26

The instantaneous fields reveal the structure of the turbulent flame [52], namely, a partially premixed flame positioned at the LOL distance where ignition kernels appear and reactants burn into products reaching temperature values close to 2000 K on the stoichiometric level curve. This premixed flame 495

is radially displaced shaping the flame as two reactive lobes separated by an inert or less reactive mass on the axis due to the richness of its mixture. As the flow moves downstream of the LOL, regions located around the stoichiometric level curve show higher temperatures due to the completeness of the chemical reactions and a diffusion flame is established. At the transient head of the spray,

500

pockets with very high and low temperatures are seen to alternate due to the inhomogeneities of the turbulent flow. Formaldehyde field shows that this species is produced in the vicinity downstream the LOL as it is a tracer of the low and intermediate temperature reactions. In addition, it is observed how this field is enclosed by the stoichiometric

505

level curve and, hence, appears in the rich mixtures region, positioned close to the axis, surrounded and partially overlapped by the partially premixed flame. A high variability is found for this field as is observed for the different instants due to the high velocities found in this region. Regarding hydroxide, it is found on the stoichiometric surface where very

510

high temperatures are reached and at the transient head of the spray where a wide spatial region shows very high temperatures. This picture of the flame structure agrees with classical works [52] and recent experimental findings [53]. Moreover, the LES simulation reproduces the small eddies found in the pe-

515

riphery of the spray due to the shear stresses, produced by the interaction between high speed spray and quiescent air, responsible of the intense mixing process [75]. A high fluctuation of these eddies is found close to the LOL and its proximities. When moving downstream, the exchange of momentum between the spray and the air due to air entrainment decreases the velocities and

520

increases the characteristic time scales of the eddies. This point is important when obtaining averaged fields and will be further commented in next section. 27

In line with the fluctuations, a high variability in the stoichiometric level curves is observed between the different instants. The fact that the hydroxide is positioned in the close vicinity of the stoichiometric level curve and the high 525

variability of this curve, induced by the turbulent motions, makes difficult to obtain a representative averaged field for Yg OH as was mentioned in previous

paragraphs. This fact may influence LOL values when computed from averaged fields if the temporal window does not extend during long times.

More insight on the structure of the flame is gained when showing the instan530

taneous fields in Z-T maps for the nominal case. These maps are represented in figure 12 for the same instants than in figure 11. For a better visualization only points with a mass fraction higher than the 0.25 of the maximum of the g field are included. Figure 12 shows the mass fraction fields for Y^ CH2 O , YOH and

Y^ C2 H2 , where acetylene (C2 H2 ) is a soot precursor.

28

2400

CH2O C2H2 OH

2200 2000

T [K]

1800 1600 1400 1200 1000 800 600 0

0.02

0.04

0.06

0.08

0.1

0.12

0.14

0.16

Z [−] 2400

CH2O C2H2 OH

2200 2000

T [K]

1800 1600 1400 1200 1000 800 600 0

0.02

0.04

0.06

0.08

0.1

0.12

0.14

0.16

Z [−] 2400

CH2O C2H2 OH

2200 2000

T [K]

1800 1600 1400 1200 1000 800 600 0

0.02

0.04

0.06

0.08

0.1

0.12

0.14

0.16

Z [−]

g Figure 12: Z-T maps for 3500 (top), 3800 (center) and 4000 µs (bottom) for the Y^ CH2 O , YOH ^ and YC2 H2 fields. The adiabatic initial mixing and upper contour of the map are included with black dashed lines. The stoichiometric value is shown with a vertical line. Results for nominal case.

535

The maps show that formaldehyde is positioned in the intermediate temperatures and rich mixtures regions while hydroxide appears in the vicinity of the stoichiometric position for very high temperatures and close to the equilibrium as was pointed in previous paragraphs. Acetylene is limited to rich mixtures and temperatures close to equilibrium.

29

540

In addition, there exists a clear variability between the position and density of the points of the fields in the maps. This is especially marked for formaldehyde as was stated when describing the spatial fields, probably due to the high turbulence levels found in the vicinity of the LOL. The upper contour of the maps shows negligible fluctuations showing that there always exist points at

545

equilibrium in the whole range of reactive mixtures. A comparison in the Z-T space, similar to that shown in figure 12, for the 750, 800 and 900 K cases is shown in figure 13 for the instantaneous fields. It is observed how the increase of reactivity by means of the ambient temperature enhances the rich mixtures reactivity and rises the maximum temperature

550

reached in the domain. The displacement of combustion towards richer mixtures is of paramount importance for soot formation. Experimentally, it has been found that negligible soot concentrations are detected for the 750 K case [3]. In the simulations it is seen that C2 H2 field spreads on richer mixtures that reach higher temperatures when increasing the ambient temperature and,

555

hence, it is expected that this would lead to higher soot formation.

30

2400

CH2O C2H2 OH

2200 2000

T [K]

1800 1600 1400 1200 1000 800 600 0

0.02

0.04

0.06

0.08

0.1

0.12

0.14

0.16

Z [−] 2400

CH2O C2H2 OH

2200 2000

T [K]

1800 1600 1400 1200 1000 800 600 0

0.02

0.04

0.06

0.08

0.1

0.12

0.14

0.16

Z [−] 2400

CH2O C2H2 OH

2200 2000

T [K]

1800 1600 1400 1200 1000 800 600 0

0.02

0.04

0.06

0.08

0.1

0.12

0.14

0.16

Z [−]

g Figure 13: Z-T maps for 750 (top), 800 (center) and 900 K (bottom) for the Y^ CH2 O , YOH and ^ YC2 H2 fields. The adiabatic initial mixing and upper contour of the map are included with black dashed lines. The stoichiometric value is shown with a vertical line. Results for very advanced instants.

The relative position of the species does not change substantially with the ambient temperature. However, it is observed how for the 750 K case the partially premixed combustion is shifted to mixtures close to stoichiometry and formaldehyde may be found at lean mixtures due to the high LOL value where 560

combustion starts developing.

31

This section is closed with an analysis and discussion of the stabilization mechanism of the flame and the apparition of ignition kernels in the LOL region. Experimental and modelling results evidence the role of auto-ignition as the main stabilization mechanism in diesel flames [36, 38, 17, 20] implying that 565

flame stabilization is strongly governed by chemistry. To perform such analysis figure 14 shows the isosurfaces for temperature at 1400 K colored by the mixture fraction for the reacting nominal case for two different time steps. The plane parallel to the plane OXZ is positioned at the LOL. It is observed that some hot spots appear upstream and detached from

570

the main flame (left plot of figure 14) that evolve in time until they join and merge with the base of the flame (right plot of figure 14). As it is deduced from figure 14 these detached hot spots appear indistinctly at the lean and rich sides.

Figure 14: Isosurfaces of temperature at 1400 K colored by the mixture fraction for the reacting nominal case for two different instants. The plane parallel to the plane OXZ denotes the position of the LOL.

To obtain a clearer picture of this phenomenon and to compare how it is affected when changing the boundary conditions, figure 15 shows the instanta575

neous planar temperature field for the cases of 750, 800 and 900 K in the region of the LOL for the modelled results for different instants. The figures include the instantaneous LOL value and the stoichiometric mixture fraction level curve. In addition, level curves for threshold values of ambient temperature plus 500 K

32

are depicted for reference.

Figure 15: Instantaneous temperature fields for Tamb = 750 (left), 800 (center) and 900 K (right). Vertical white line shows the instantaneous LOL position while the solid white line indicates the stoichiometric level curve. Black lines are level curves for Tamb + 500 K.

580

The figures show how the level curves for temperature separate two zones of unburned and burnt mixtures with a well-defined interface due to the sharp jump in this field as a consequence of the high chemical source terms found in the partially premixed combustion. For the three temperatures it is observed how small ignition kernels appear

585

in the vicinity of the LOL (for a better visualization different rectangles have been drawn surrounding some of these kernels). These pockets, which appear

33

isolated and upstream of the main flame, increase in size and are convected downstream until they are brought close enough to get attached to the main flame. During the whole process the pocket expands over wider regions due to 590

its growth to the surroundings as a consequence of the heat release and species diffusion. The fact that such hot kernels are incepted spontaneously upstream and detached from the base of the main flame is an evidence that auto-ignition is one of the flame stabilization mechanisms and plays an important role. In fact, taking

595

into account the position of the stoichiometric contour, the ignition kernels are detected both on the lean and rich sides. This result agrees with recent DNS simulations from Tagliante et al. [76], which indicate that autoignition kernels occur in terms of isolated autoignition or assisted autoignition, the former being related respectively to rich and lean locations, the latter one assisted by the

600

upstream entrainment of burnt gases upstream from the reaction front. This phenomenology agrees with experimental observations where similar detached reactive pockets were identified and led to conclude that the stabilization mechanism in a diesel flame is basically controlled by auto-ignition [36]. In addition, modelling works [38, 17, 34, 20, 76] have reported similar observations

605

about the formation of ignition kernels, their expansion and attachment to the main flame. In addition, it is worth mentioning how pockets of unburned mixture are observed inside the flame as a consequence of its motion. It is observed how two branches of high temperature region get closer until they eventually merge and

610

enclose a volume of fresh mixture. This is seen in the last instants for the cases 750 and 800 K in figure 15 where again for a better visualization these pockets have been enclosed in rectangles. Moreover, figure 15 shows clearly how the morphology of the base of the flame is drastically modified when changing the ambient temperature. Different

615

to the nominal case, where two reactive lobes are clearly identified and the center of the flame remains inert, the 750 K case shows that the length of the lobes decrease while their width increase. Hence, reducing the ambient temperature 34

provokes that the base of the flame becomes flatter. Finally, this analysis is concluded showing the flame index [77] which is 620

defined by the expression F I = ∇Yef · ∇Yeox

(14)

where Yef and Yeox are the instantaneous fuel, C12 H26 , and oxidant, O2 , mass

fractions and provides a measure, depending on its sign, of the predominance

of the premixed (F I > 0) or non-premixed (F I < 0) combustion in the flame. Figure 16 shows a meridian plane for the flame index for the nominal case for 625

an advanced time step. −3

5

x 10

200

r [m]

100 0

0 −100

−5 0.01

−200 0.015

0.02 x [m]

0.025

0.03

Figure 16: Flame index for Tamb = 900 K at an advanced time step. The solid white line corresponds to the stoichiometric mixture fraction level curve and the dashed white line indicates the position of the LOL.

In both cases it is observed that the flame is composed of a region where non-premixed combustion predominates which extends over the vicinity of the stoichiometric mixture fraction, the whole region upstream the LOL and some inner regions downstream the LOL. On the contrary, the region where premixed 630

flames are more important is concentrated on the inner region of the flame downstream the LOL. This means that in a diesel flame, auto-ignition, found upstream the LOL, occurs in the non-premixed combustion mode in agreement with DNS results [78, 79] and other works devoted to spray A [34].

35

4.2.3. Averaged fields 635

The analysis of the reactive spray is closed with a comparison of the averaged fields for the ambient temperature sweep for very advanced instants. As the main aspects of the flame have been already described in the previous sections the spray behaviour is here briefly summarised. Some comments about the averaging process are first explained.

640

One of the major difficulties in averaging the fields of diesel sprays is that even for the long time interval for which the simulation extends (4 ms) there exists a part of the head of the spray that is not in quasi-steady conditions. In addition, some fields like hydroxide found close to the stoichiometric level curve show strong variations complicating the averaging process.

645

In the current work, fields have been averaged in both azimuthal and temporal directions to obtain more representative averages. After checking, it was confirmed that 32 meridian planes are enough to reach convergence for the azimuthal averages. As the azimuthal average has slight effect on regions close to the axis and no effect on the axis, a temporal average is considered too. The

650

window for averaging has been selected between 3.5 and 4 ms in order to consider a developed spray that extends in the wider possible region taking into account the limitations of the simulation. As the spray shows a quasi-steady evolution [53], i.e. the spray fluid dynamics in a region stabilize rapidly once the spray passes through it, it is expected

655

that representative averages are obtained in almost the whole spray and only the transient evolution affects the most advanced regions of the jet. Clearly, the representativeness of the average depends on the own nature of the field too, i.e. how the field changes with the mixture fraction, the scalar dissipation rate, etc.

660

The spatial averaged fields for the cases with ambient temperature 750, 800 g and 900 K are shown in figure 17 for hTei, hY^ CH2 O i and hYOH i fields.

36

g Figure 17: Averaged fields hTei (first row), hY^ CH2 O i (second row) and hYOH i (third row) for 750 (left), 800 (center) and 900 K (right) cases. The averaged fields correspond to developed e = Zst is included in white line. LOL values for experiment (white sprays. The level curve hZi dashed) and simulation (red dashed) are shown too.

From figure 17 it is observed that the ambient temperature has a strong impact on the flame structure as expected. Reducing the ambient temperature moves LOL downstream displacing, hence, combustion to leaner mixtures and 665

forcing that the whole combustion takes place in more reduced spatial regions since the stoichiometric level curve is closed at approximately the same distance from the nozzle for all the considered cases. The expected flame shape composed by two lobes is observed in all the cases although, as previously mentioned, this is not so clearly seen when reducing the ambient temperature since the

670

displacement of combustion downstream allows that the mixtures on the axis react at the same distance that the periphery of the spray. As previously described, high concentrations of formaldehyde are observed in the vicinity of LOL and enclosed in the region defined by the stoichiometric level curve, since this region corresponds to intermediate temperatures and rich

675

mixtures. On the contrary, hydroxide is positioned on the stoichiometric level curve spreading in the slightly lean region where very high temperature values are reached. These averaged fields, which seem to converge for short and intermediate distances (until 60 mm), agree with those obtained by means of RANS simula-

680

tions in other works [34, 10]. However, for larger distances, around 80 mm, it is probable that they have not converged completely as seen in the hydroxide

37

field and the stoichiometric level curve.

5. Conclusions and future work Spray A from ECN has been modelled by means of LES simulations in in685

ert and reactive conditions for an ambient temperature sweep. The Dynamic Structure model has been used to model turbulence while a simplified approach of the flamelet concept has been applied for the turbulent combustion evolution. A well-known chemical mechanism for the dodecane oxidation has been considered.

690

The inert spray has been calibrated in terms of vapour penetrations and mixture fraction and velocity fields. In addition, a satisfactory agreement has been obtained for the mixture fraction fluctuations. The spectral decomposition of the auto-correlation functions show that the LES simulations reproduce the inertial sub-range. In addition, the tip penetration and the velocity field for

695

reactive conditions are properly described by the simulation. Finally, qualitative agreement with experiments is observed for CH2 O and OH fields for the reactive flame. Regarding the reactive simulations, it was observed that the ID is wellcaptured except for very low ambient temperatures while the LOL trends are

700

correctly described although the absolute value is under-predicted. In addition, LOL has been shown to be subjected to a weak variability for all the cases, pointing in the direction that chemistry is intense enough to absorb the local dynamic flow fluctuations. The instantaneous fields reproduce a partially premixed flame structure

705

shaped in two lobes in the vicinity of the LOL followed by a diffusion flame established downstream where high temperature pockets are alternated by regions with lower temperatures as a consequence of the inhomogeneities that the turbulent flow induces. Formaldehyde is found in a close region to the LOL and the axis, showing some variability in Z-T maps, while hydroxide is positioned

710

on the stoichiometric surface level and the lean side. The ambient temperature

38

has a strong effect on the position of the flame and, consequently, on the range of reactive mixtures with a clear impact on the soot formation. For the different boundary conditions it has been found that some ignition kernels appear in the proximity of the LOL, grow and merge with the main flame. This is an evidence 715

of auto-ignition as a flame stabilization mechanism which is supported by the negative values of the flame index reached in the vicinity of the LOL. Averages have been carried out in azimuthal and temporal directions in order to filter the fluctuations related to the flow. An important fraction of the spray was correctly averaged because it reaches a quasi-steady regime and only at the

720

head of the transient spray it was observed that the temporal interval was too short to provide representative averages. Moreover species with high variability appearing in the high temperature region, such as hydroxide, were especially difficult to average. As a final conclusion, this work shows that the LES simulations provide rea-

725

sonable results and the flamelet concept and, in particular, the ADF approach, has the ability to describe the flame structure accurately. This investigation encourages to analyse the influence of other parameters, such as the oxygen concentration, and to study in detail the stabilization flame mechanism.

6. Acknowledgements 730

Authors acknowledge that this work was possible thanks to the Ayuda para la Formaci´ on de Profesorado Universitario (FPU 14/03278) belonging to the Subprogramas de Formaci´ on y de Movilidad del Ministerio de Educaci´ on, Cultura y Deporte from Spain. Also this study was partially funded by the Ministerio de Econom´ıa y Competitividad from Spain in the frame of the COMEFF

735

(TRA2014-59483-R) national project. Finally, the authors thankfully acknowledge the computer resources at MareNostrum and technical support provided by Barcelona Supercomputing Center (RES-FI-2017-2-0044).

39

References [1] R. Novella, A. Garc´ıa, J. Pastor, V. Dom´enech, Mathematical and Com740

puter Modelling 54 (2011) 1706–19. [2] L. M. Pickett, G. Bruneaux, R. Payri, Sandia National Laboratories, Livermore, CA, https://ecn.sandia.gov/ (2019). [3] N. Maes, M. Meijer, N. Dam, B. Somers, H. B. Toda, G. Bruneaux, S. Skeen, L. Pickett, J. Manin, Combustion and Flame 174 (2016) 138–

745

51. [4] M. Bardi, R. Payri, L. Malbec, G. Bruneaux, L. Pickett, J. Manin, T. Bazyn, C. L. Genzale, Atomization and Sprays 22 (2012). [5] J. Benajes, R. Payri, M. Bardi, P. Mart´ı-Aldarav´ı, Applied Thermal Engineering 58 (2013) 554–63.

750

[6] R. Payri, J. Garc´ıa-Oliver, T. Xuan, M. Bardi, Applied Thermal Engineering 90 (2015) 619–29. [7] L. Cutrone, P. De Palma, G. Pascazio, M. Napolitano, Computers & Fluids 39 (2010) 485–98. [8] B. Naud, R. Novella, J. Pastor, J. Winklinger, Combustion and Flame 162

755

(2015) 893–906. [9] Y. Pei, E. Hawkes, S. Kook, G. Goldin, T. Lu, Combustion and Flame 162 (2015) 2006–19. [10] J. Desantes, J. Garc´ıa-Oliver, R. Novella, E. P´erez-S´ anchez, Applied Thermal Engineering 117 (2017) 50–64.

760

[11] A. Kolmogorov 30 (1941) 299–303. [12] S. Pope, Turbulent flows, IOP Publishing, 2001. [13] S. Pope, New journal of Physics 6 (2004) 35. 40

[14] H. Pitsch, Annu. Rev. Fluid Mech. 38 (2006) 453–82. [15] T. Poinsot, D. Veynante, Theoretical and numerical combustion, RT Ed765

wards Inc., 2005. [16] K. J. Range Dinesh, A. Savill, K. Jenkins, M. Kirkpatrick, Computers & Fluids 39 (2010) 1685–95. [17] Y. Pei, S. Som, E. Pomraning, P. Senecal, S. Skeen, J. Manin, L. Pickett, Combustion and Flame 162 (2015) 4442–55.

770

[18] A. N. Toosi, M. Farokhi, B. Mashadi, Computers & Fluids 109 (2015) 85–99. [19] S. Buhl, F. Dietzsch, C. Buhl, C. Hasse, Computers & Fluids 156 (2017) 66–80. [20] H. Kahila, A. Wehrfritz, O. Kaario, M. Masouleh, N. Maes, B. Somers,

775

V. Vuorinen, Combustion and Flame 191 (2018) 142–59. [21] M. Germano, U. Piomelli, P. Moin, W. Cabot, Physics of Fluids A: Fluid Dynamics 3 (1991) 1760–5. [22] A. Yoshizawa, K. Horiuti, Journal of the Physical Society of Japan 54 (1985) 2834–9.

780

[23] S. Ketterl, M. Klein, Computers & Fluids 165 (2018) 64–77. [24] E. Pomraning, C. Rutland, AIAA journal 40 (2002) 689–701. [25] N. Bharadwaj, C. Rutland, S. Chang, International Journal of Engine Research 10 (2009) 97–119. [26] J. Momp´ o Laborda, Engineering Large Eddy Simulation of Diesel Sprays,

785

Ph.D. thesis, Universitat Polit`ecnica de Val`encia, 2014. [27] T. Lucchini, G. d’Errico, D. Ettorre, G. Ferrari, SAE International Journal of Fuels and Lubricants 2 (2009) 966–75. 41

[28] E. Fooladgar, C. K. Chan, K. J. Nogenmyr, Computers & Fluids 146 (2017) 42–50. 790

[29] S. Bhattacharjee, D. Haworth, Combustion and Flame 160 (2013) 2083– 102. [30] H. Barths, C. Hasse, G. Bikas, N. Peters, Proceedings of the Combustion Institute 28 (2000) 1161–8. [31] G. D’Errico, T. Lucchini, F. Contino, M. Jangi, X. S. Bai, Combustion

795

Theory and Modelling 18 (2014) 65–88. [32] P. Kundu, Y. Pei, M. Wang, R. Mandhapati, S. Som, Atomization and Sprays 24 (2014). [33] A. Wehrfritz, O. Kaario, V. Vuorinen, B. Somers, Combustion and Flame 167 (2016) 113–31.

800

[34] Y. Pei, E. R. Hawkes, M. Bolla, S. Kook, G. M. Goldin, Y. Yang, S. Pope, S. Som, Combustion and Flame 168 (2016) 420–35. [35] C. A. Idicheria, L. M. Pickett, SAE Technical Paper (2006). [36] L. Pickett, D. Siebers, C. Idicheria, SAE technical paper (2005). [37] F. Tagliante, L. M. Malbec, G. Bruneaux, L. M. Pickett, C. Angelberger,

805

Combustion and Flame 197 (2018) 215–26. [38] C. Gong, M. Jangi, X. Bai, Applied Energy 136 (2014) 373–81. [39] J. v. Oijen, L. d. Goey, Combustion Science and Technology 161 (2000) 113–37. [40] C. Pierce, P. Moin, Journal of fluid Mechanics 504 (2004) 73–97.

810

[41] M. Ihme, H. Pitsch, Combustion and Flame 155 (2008) 70–89. [42] M. Ihme, H. Pitsch, Combustion and Flame 155 (2008) 90–107.

42

[43] J. Tillou, J. Michel, C. Angelberger, C. Bekdemir, D. Veynante, Oil & Gas Science and Technology-Revue de l’IFP 69 (2014) 155–65. [44] F. Williams, Recent advances in theoretical descriptions of turbulent diffu815

sion flames, Springer, 1975. [45] N. Peters, Progress in Energy and Combustion Science 10 (1984) 319–39. [46] N. Peters, Turbulent combustion, Cambridge University Press, 2000. [47] O. Gicquel, N. Darabiha, D. Th´evenin, Proceedings of the Combustion Institute 28 (2000) 1901–8.

820

[48] U. Maas, S. Pope, Combustion and Flame 88 (1992) 239–64. [49] J. Michel, O. Colin, D. Veynante, Combustion and Flame 152 (2008) 80–99. [50] J. Michel, O. Colin, International Journal of Engine Research 15 (2014) 346–69. [51] J. Tillou, J. Michel, C. Angelberger, D. Veynante, Combustion and Flame

825

161 (2014) 525–40. [52] J. Dec, SAE Technical paper (1997). [53] J. Garc´ıa-Oliver, L. Malbec, H. Toda, G. Bruneaux, Combustion and Flame 179 (2017) 157–71. [54] A. Kastengren, F. Tilocco, C. Powell, J. Manin, L. Pickett, R. Payri,

830

T. Bazyn, Atomization and Sprays 22 (2012) 1011–52. [55] Universitat

Polit`ecnica

de

Val`encia,

Spain,

http://www.cmt.upv.es/ECN03.aspx (2019). [56] J. Pastor, J. Garcia-Oliver, J. Pastor, W. Vera-Tudela, Atomization and Sprays 25 (2015). 835

[57] http://www.openfoam.com/ (2019).

43

[58] P. Senecal, E. Pomraning, K. Richards, S. Som, SAE Technical Paper (2013). [59] Q. Xue, S. Som, P. Senecal, E. Pomraning (2013). [60] N. Bharadwaj, C. Rutland, Atomization and Sprays 20 (2010). 840

[61] J. Beale, R. Reitz, Atomization and Sprays 9 (1999). [62] P. Nordin, Complex chemistry modeling of diesel spray combustion, volume 18, Chalmers University of Technology, Sweden, 2001. [63] http://www.reactiondesign.com/products/chemkin/ (2019). [64] F. Payri, R. Novella, J. Pastor, E. P´erez-S´ anchez, Applied Mathematical

845

Modelling 49 (2017) 354–74. [65] Y. Chen, M. Ihme, Combustion and Flame 160 (2013) 2896–910. [66] C. Pera, O. Colin, S. Jay, Oil & Gas Science and Technology-Revue de l’IFP 64 (2009) 243–58. [67] K. Narayanaswamy, P. Pepiot, H. Pitsch, Combustion and Flame 161

850

(2014) 866–84. [68] A. Frassoldati, G. D’Errico, T. Lucchini, A. Stagni, A. Cuoci, T. Faravelli, A. Onorati, E. Ranzi, Combustion and Flame 162 (2015) 3991–4007. [69] B. Akkurt, H. Akargun, L. Somers, N. Deen, R. Novella, E. P´erez-S´ anchez, SAE Technical Paper 2017-01-0577, 2017, doi:10.4271/2017-01-0577 (2017).

855

[70] L. Pickett, J. Manin, C. Genzale, D. Siebers, M. Musculus, C. Idicheria, SAE International Journal of Engines 4 (2011) 764–99. [71] M. Meijer, L. Malbec, G. Bruneaux, L. Somers, in: 12th Triennial International Conference on Liquid Atomization and Spray Systems (ICLASS 2012), Heidelberg, Germany, September, pp. 2–6.

44

860

[72] W. Eagle, M. Musculus, L. Malbec, G. Bruneaux, in: ILASS Americas 26th Annual Conference on Liquid Atomization and Spray Systems, Portland, OR. [73] H. Tennekes, J. Lumley, A first course in turbulence, MIT press, 1972. [74] S. Subramaniam, Progress in Energy and Combustion Science 39 (2013)

865

215–45. [75] C. Duwig, L. Fuchs, Combustion Science and Technology 180 (2008) 453– 80. [76] F. Tagliante, T. Poinsot, L. M. Pickett, P. Pepiot, L. M. Malbec, G. Bruneaux, C. Angelberger, Combustion and Flame 201 (2019) 65–77.

870

[77] H. Yamashita, M. Shimada, T. Takeno 26 (1996) 27–34. [78] P. Domingo, L. Vervisch, J. R´eveillon, Combustion and Flame 140 (2005) 172–95. [79] C. S. Yoo, R. Sankaran, J. H. Chen, Journal of Fluid Mechanics 640 (2009) 453–81.

45