Impact of numerical method on auto-ignition in a temporally evolving mixing layer at various initial conditions

Numerical analysis of the auto-ignition of turbulent mixing layer between the cold fuel (hydrogen) and hot oxidizer (air) is presented. The research were performed using an Implicit-Large Eddy Simulation (ILES) method with attention on auto-ignition time, flame kernel localisation and propagation. We focused on an impact of discretization method on auto-ignition scenario and flame development. The results obtained showed that numerical approach plays an important role and to some extent may falsify the results, especially for low oxidiser temperatures.


Introduction
Turbulent combustion process, being a result of the two-way interaction of chemical reactions and turbulent structures, is characterized by significant oscillations of all flow variables. A flame interacting with turbulent flow leads to its modification caused by strong accelerations through the flame front and because of large changes of density and viscosity associated with temperature variations. An opposite impact is equally important, i.e. the turbulence affects the combustion process by enhancing the chemical reactions or inhibiting them, which may lead to a spontaneous ignition or extinction. The fact that accuracy of modelling of turbulent flows is directly related to an applied model and numerical method indicates that simulations of reacting flow are very complex problem and are regarded as one of the most difficult tasks in the contemporary CFD (Computational Fluid Dynamics). They need an involvement of computational methods capable to accurately represent a wide range of strongly unsteady phenomena occurring in the combustion process.
During the last decades numerous works devoted to DNS (Direct Numerical Simulation) of autoignition were presented. One of the first numerical analysis of auto-ignition (DNS calculations with one-step chemistry) was performed for two dimensional (2D) mixing layer (methane/air) by Mastorakos et al. [1] in 90s. This fundamental work constituted the basis for subsequent research in the field of numerical analyses of the auto-ignition phenomenon. Recently, the auto-ignition in 3D case was analysed using DNS by Doom and Mahesh [2]. They studied a configuration similar as in [1] and found that extending the mixing layer over the third dimension three possible auto-ignition scenarios can be observed, not reported for 2D cases.
The present work concerns Implicit Large Eddy Simulation (ILES) of non-premixed combustion in a temporally evolving mixing layer between cold fuel (hydrogen) and hot oxidiser (air). Similarly to [2], on the initial velocity field we impose an isotropic turbulence. The analysis is focused on the auto-ignition and flame propagation phenomena and their dependence on discretization method, which in ILES approach plays a crucial role.

Governing equations
We consider low Mach number reactive flows governed by the Navier-Stokes equations, the transport equations for each of reactive species and energy equations complemented with the equation of state. In the context of LES method the set of the governing equations is given as: The bar and tilde symbols denote filtered quantities [3].
and for a grid resolution sufficiently close to the DNS requirements the above assumption can be considered as appropriate.

Solution algorithm
The computations are carried out using an in-house numerical code (SAILOR) based on the low Mach number approach which was previously well verified in gaseous shear layer flows with/without chemical reactions (see [6][7][8][9]). The spatial discretization is performed on half-staggered meshes by high-order (6 th ) compact difference approximation for the Navier-Stokes and continuity equations while for the chemical transport equations two different schemes are examined: 2 nd order TVD (Total Variation Diminishing) with van Leer's limiters and the 5 th order WENO (Weighted Essentially Non-Oscillatory). The governing equations are integrated in time using the 2 nd order predictor-corrector approach combined with the projection method for pressure-velocity coupling. In the half-staggered approach the velocity components and scalar quantities are stored in the same spatial locations (nodes) while the pressure nodes are shifted to the cell centres. Compared to the fully staggered grid arrangement, where the pressure and all scalar variables are shifted to the cell centres and the velocity components are placed on different cell corners, the use of half-staggered meshes facilitates the solution algorithm as the velocity needs to be interpolated only within the projection method [6,7].

Simulation details
The calculations are performed to study the auto-ignition in the temporally evolving mixing layer between two separate streams (cold fuel/hot air) moving with opposite velocities. The initial isotropic turbulence prescribed by a three dimensional kinetic energy spectrum [2]: was combined with the Gaussian function and superimposed on the initial flow field providing the largest velocity fluctuations in a region of fuel/oxidizer contact zone (see figure 1). In the equation (2), is the RMS of initial velocity fluctuations and 0 k is used to adjust the turbulent length scale . We consider the cases with the turbulence intensity

Results and discussion
In all analyzed cases the simulations started from the initial conditions and continued until the autoignition occurred and the flame started to propagate. The obtained results clearly reveal influence of discretization schemes on the auto-ignition process and flame propagation phase. The results presented in the figures 3-9 show the time evolutions of maximum quantities in the domain (HO 2 mass fraction and temperature), 2D contours of instantaneous temperature distributions as well as 3D contours with iso-surfaces inside the flame. Figures 5-8 show the enlarged regions of the shear layer at selected time moments.
During the calculations the maximum values of temperature and species mass fractions were registered and their sudden increase indicated the auto-ignition. Figures 3 and 4 (figure 3). The profiles of mass fraction of OH radical, not presented here, starts to grow at the time instant close to the auto-ignition time indicated based on the temperature growth, which is consistent with literature [11] studies of autoignition processes. On the other hand the radical HO 2 , the so-called pre-ignition species, appears soon before the temperature growth, which is well seen in both the figures. The results presented show that for the higher oxidiser temperature the profiles slightly depend on the discretization method and for the TVD scheme (dashed lines) their time evolution is smoother. Comparing the max T profiles, the main differences are observed shortly after the auto-ignition when the flame starts to grow. The evolution of maximum HO 2 mass fraction is similar at the beginning of auto-ignition, however the profiles start to diverge at t = 1.35 ms. Figure 4 shows that in the case when the oxidizer temperature is reduced to 970K, a significant impact of numerical approach becomes more pronounced both during the auto-ignition and propagation phases. The maximum values of HO 2 and temperature obtained using the WENO scheme grow faster, achieve very high level and reveal some oscillations, contrary to rather smooth profiles obtained using TVD scheme. Similar observations can be made analyzing the evolution of temperature in enlarged regions of 2D cross-sections located in the centre of 3D computational domain, shown in the figure 5. The time instances selected correspond to the subsequent flame development phases. It is seen that the flame develops starting from many separated ignition kernels distributed along the mixing layer. At the time instant t = 1.8 ms, corresponding to the rise of the max T profiles in the figure 4, the differences observed in the temperature distributions, for both schemes are evident. Analysis of the contour plots reveals that the numerics affects both the appearance, distribution and evolution of the flame kernels. It seems that for lower oxidizer temperature, which corresponds to longer simulation time needed to the auto-ignition occurrence, which in turn enhances the mixing of fuel and oxidizer, the convective and diffusive mechanisms play an important role.   Comparing the figures 7-8 one may conclude that the main differences between the obtained solutions using TVD and WENO schemes lie in the smoother resolution of the flow field obtained using the TVD approach. Concerning the results for the flow conditions at 100 Re   , the flame kernels localization and maximum values of variables describing the flow field seem to be also significantly dependent on the applied discretization scheme. In the figure 6 it is seen that the flame develops from single sites for TVD, whereas when the WENO is used, the flame appears in a large region of the mixing layer. The profiles of maximum temperature and HO 2 mass fraction, presented in the figure 9, confirm assumptions that the use of various numerical approaches may lead to large discrepancies and falsification of the results. Future studies will include a mesh sensitivity analysis. One could expect that when the ILES mesh resolution approaches a mesh resolution required for DNS an impact of the numerical method becomes definitively smaller.

Conclusions
The paper presented the numerical analysis of temporally evolving mixing layer with the particular attention on the numerical impact of discretization method on auto-ignition time, localization of ignition kernels and flame propagation. The results were obtained using the TVD and WENO schemes. We analyzed the mixing layer between two separate streams (cold hydrogen/hot air) moving in opposite directions. As a initial fluctuation an isotropic turbulence was prescribed in the region of the shear layer. It was found that in the cases when the oxidizer temperature was low (970K) the numerical approach significantly affected the ignition scenario and growth of the flame. On the other hand, when the oxidiser temperature was higher (980K) the numerical method did not play such an important role and it seems that in this case the chemical kinetics was a crucial element of the simulation. However, this observations are valid only for specified initial conditions.