Control of a Lifted H2/N2 Flame by Axial and Flapping Forcing: A Numerical Study


 The large eddy simulation (LES) method combined with the Eulerian stochastic field approach has been used to study excited lifted hydrogen flames in a stream of hot co-flow air in a configuration closely corresponding to the so-called Cabra flame. The excitation is obtained by adding to an inlet velocity profile three types of forcing ((i) axial; (ii) flapping; (iii) combination of both) with amplitude of 15% of the fuel jet velocity and frequency corresponding to the Strouhal numbers St=0.30, 0.45, 0.60 and 0.75. It is shown that such a type of forcing significantly changes the lift-off height Lh of the flame and its global shape, resulting in a flame occupying large volume or the flame, which downstream the nozzle transforms from the circular one into a quasi-planar flame. Both the Lh and their spreading angles of the flames were found to be a function of the type of the forcing and its frequency. The minimum value of Lh has been found for the case with the combination of axial and flapping forcing at the frequency close to the preferred one in the unexcited configuration. The impact of the flapping forcing manifested through a widening of the flame in the flapping direction. It was shown that the excitation can significantly increase the level of the velocity and temperature fluctuations intensifying the mixing process. The computational results are validated based on the solutions obtained for a non-excited flame for which experimental data are available.


Introduction
Intensive research on applications of active flow control methods was initiated by a famous work of Crow and Champagne 1 concerning a round jet. Applying a low amplitude excitation (forcing) being a sinusoidal function of time they observed that for some range of frequency the jet behaviour significantly changes. The turbulence intensity level expressed in terms of the velocity fluctuations increased, and moreover, its profile along the jet axis was characterised by a distinct local maximum at a distance of approximately eight nozzle diameter from the nozzle exit. Such a behaviour was never reported before and did not occur in natural jets. It turned out that applying a relatively simple control technique one may change the flow dynamics to the extent incomparably larger to an energy input needed to introduce the excitations. These findings stimulated very extensive analyzes of various excitation types including axial, flapping and helical forcing modes [2][3][4] , and their combinations. Recent LES (Large Eddy Simulations) and DNS (Direct Numerical Simulations) works of Tyliszczak A. 5,6 and Gohil et al. 7 showed that for carefully selected frequencies of combined axial and helical excitation one can obtain multi-armed jets, i.e., the jets characterized by 5, 7 or 11 or even 20 separate branches, closely reminding the blooming jets reported by Reynolds and Parekh 8 . Without doubts, the active flow control methods are superior compared to the passive methods relying on optimization of a flow domain for particular flow regimes. Mainly, because their driving parameters (the excitation amplitude/frequency, spatial distribution, etc.) can be dynamically adapted to changing flow conditions, e.g., increasing/decreasing inlet velocity or temperature. The passive methods do not allow for such freedom, yet, from the point of view of the working costs, they are certainly cheaper as they do not need additional energy to operate.
The findings on the active control techniques quickly translated to combustion science where they have been the focus of interest since the early 1990s. 9,10 . Regarding the fundamental problems the attention is very often paid to jet-type flames in which the mechanical or acoustic forcing acts as a source of external excitation. The former is usually introduced by specially designed nozzle tips with magnetic or piezoelectric actuators 11,12 . The acoustic excitation is more often used and is added by loudspeakers mounted upstream of the nozzle exits [13][14][15][16][17] . The influence of this type of excitation on reduction of pollution emissions in a lean premixed lifted flames and flame stability was demonstrated by Chao et al. 13,14 , among others. Focussing on stability issues, they found that the excitation significantly alters the flame dynamics and can be used as a "tool" suppressing or amplifying the stabilization process. Abdurakipov et al. 15 demonstrated that in comparison to natural flames the excitation ensures stable combustion regimes and visibly shifts the blow-off limits. A research on an excited lifted non-premixed flame in a hysteresis regime, i.e., when depending on initial conditions a flame can be attached to a nozzle or remain lifted for the same fuel velocity, was performed by Demare and Baillot 16,17 . It was shown that, by changing the amplitude and frequency of excitation one can enhance the combustion process or produce large fluctuations, and thus, weaken the flame stability. Kozlov et al. 18 analyzed micro-flames in the field of transverse sound waves. They found that the excitation can flatten the round flames and transform them to nearly plane flames. Surprisingly, for particular forcing frequencies the excitation led to a splitting of the flame into two separate branches in a very similar way as observed for bifurcating non-reacting, constant density jets 8 . More recently, the occurrence of the bifurcating phenomenon in flames was reported in numerical studies of Tyliszczak 19 focused on a low Reynolds number hydrogen flame. Application of only the flapping excitation caused that the flame changed its initial circular shape into the planar one with two co-existing separate arms. Moreover, for some range of the excitation frequency, a triple-flame occurred.
As discussed above, the excitation can alter the flame dynamics and influence the pollution emission. In the present work we focus on the global impact of the excitation on the flame by applying three different excitation types: (i) axial (ii) flapping; (iii) axial plus flapping with different forcing frequencies. The basic flow configuration closely corresponds to the so-called Cabra flame 20 at the Reynolds number equal to 23600. We apply LES method in combination with the Eulerian Stochastic Field (ESF) approach 21 for the combustion modelling. There are no experimental results for the excited Cabra flame, and hence, the credibility of the obtained results is proven by comparison with the measurements data available for the non-excited flame. The present research is an exploratory numerical study in which we assess large-scale effects of the excitation, i.e., the change of the lift-off height or the change of the size and shape of the flame.

Mathematical approach
We consider a low Mach number reacting flow for which the continuity, Navier-Stokes and transport equations of scalars within the framework of the LES method are defined as: where the bar and tilde symbols denote filtered quantities. The variables in Eqs. (1) are the velocity vector u, the density ρ, the hydrodynamic pressure p, the species mass fractions Y α and enthalpy h. The subscript α represents the index of the species α = 1, . . . , N-species. The quantities τ and D α , D are the viscous stress tensor and mass and heat diffusivities. The sub-filter tensor is given by τ SGS = 2µ t S, where S is the rate of strain tensor of the resolved velocity field and µ t is the sub-filter viscosity modelled as in 22 . The sub-filter diffusivities in the species and enthalpy transport equations are computed as D SGS = µ t /(ρσ ) where σ is the turbulent Schmidt or Prandtl number assumed equal to 0.7 23 . The set of equations (1) is complemented with the equation of state p 0 = ρR T with p 0 being the constant thermodynamic pressure and R is the gas constant.
The chemical source terms ρẇ α represent the net rate of formation and consumption of the species. A highly non-linear nature of this term means that sub-grid fluctuations cannot be ignored. In the present work the scalar equations (species and enthalpy) are replaced by an equivalent evolution equation for the density weighted filtered PDF function, which is solved using the stochastic field method proposed by Valiño 21 . Each scalarφ α is represented by 1 ≤ n ≤ N s stochastic fields ξ n α such The stochastic fields evolve according to: where the total diffusion coefficient is defined as Γ = D α + D SGS α , the micro-mixing time scale equals to τ SGS =ρ∆ 2 /(µ + µ t ) cell being the LES filter width and dW represents a vector of Wiener process increments different for each field. Following Jones and Navarro 24 eight stochastic fields have been used. The test computations performed with sixteen fields did not show any substantial changes in the flames dynamics.

Numerical method
We apply an in-house numerical code (SAILOR) based on the projection method for pressure-velocity coupling 25 . The time integration is performed by means of an operator splitting approach where the transport in physical space and chemical terms are solved separately. The convective and diffusive parts of the governing equations are advanced in time using a predictor-corrector technique with the 2nd order Adams-Bashforth / Adams-Moulton methods. The chemical reactions are computed using CHEMKIN interpreter. In the present study analyze the hydrogen combustion using a detailed mechanism of Mueller 26 involving 9 species and 21 reactions. The reaction terms are stiff and therefore they are integrated in time applying the VODPK 27 solver that is well suited for stiff systems. The spatial discretization is performed on half-staggered meshes applying the 6th order compact difference approximation for the Navier-Stokes and continuity equations 25,28 . The convective terms in the stochastic field equations equations are discretized applying TVD scheme with van Leer limiters. The applied code has been thoroughly verified in previous studies 19,25,[29][30][31] and it turned out to be very accurate.

Computational configuration
The basic flow configuration analyzed in this study corresponds to the so-called Cabra flame 20 , which we modify by adding the excitation at the nozzle exit. A schematic view of the computational geometry is shown in Figure 1. It is a rectangular box with dimensions L x × L y × L z = 14D × 30D × 14D, where D = 0.00457 m is the nozzle diameter. The injected fuel (X H2 = 0.254, X N2 = 0.746) has the temperature T fuel = 305 K and the bulk velocity U j = 107 m/s. Outside of the fuel nozzle there is a hot (T cf = 1045 K) co-flowing stream of the hydrogen combustion products (X O2 = 0.147, X H2O = 0.1, X N2 = 0.753) with the velocity U cf = 3.5 m/s. The excitation (forcing) is introduced as a component of the velocity prescribed at the inlet as is the mean velocity profile corresponding to the fully developed pipe flow (1/7 profile) and u turb ( x,t) = 0.05U j represents turbulent fluctuations computed applying a digital filtering method proposed by Klein et al. 32 . This method guarantees properly correlated velocity fields which reflect realistic turbulent flow conditions. The forcing component u excit ( x,t) is added to the streamwise velocity only and it is defined as: Axial forcing which is the superposition of axial and flapping forcing term with amplitudes A a and A f and frequencies f a and f f .

3/10
As mentioned in the Introduction section, the simulations and results discussed in this study has an exploratory character as no experiments or numerical data are available for validation of the obtained results regarding the excited flames. Therefore the comparison was performed only based on measurements for the original Cabra flame configuration of which a spatio-temporal complexity is not much different from the cases with the excitation. As will be presented the agreement between the preset results and experimental data is sufficiently good to assume that the obtained results are reliable.

Results and discussion
Three-dimensional flame behavior The fuel issuing from the nozzle mixes with the co-flowing hot stream and auto-ignites. The ignition spots appear at the locations of the most reactive mixture fraction ξ MR = 0.053 at distances far from the inlet, i.e. y/D ≈ 20. Then, the flame spread radially, propagates upstream and stabilizes as a lifted flame in between y/D ≈ 7.5 − 11.0 depending on the test case. Figure 2 shows iso-surfaces of instantaneous temperature (T = 1200 K) and Q-parameter (Q = 10 s −2 ) for the cases without the excitation and for A 45 , A 60 and A 75 . The Q-parameter is commonly used to indicate organized vortical motion. It is defined as Q = 1/2(Ω i j Ω i j − S i j S i j ), where S i j and Ω i j are the symmetric and antisymmetric parts of the velocity gradient tensor. Here, it visualizes the effect of the excitation, which manifests by occurrence of toroidal vortices formed in the vicinity of the nozzle. They mutually interact through rib-like vortices and break-up further downstream. Compared to the unexcited case (see Figure 2a) the vortices are much more pronounced, the distances between them are dependent on the forcing frequency and decrease with increasing St a . One can notice that in the case A 60 the jet at y/D ≈ 1.0 − 4.0 is wider than for A 45 and A 75 and its shape reminds a barrel. The temperature distributions reveal that the flame lift-off height depends on the forcing frequency and turns out to be the smallest for A 60 that will be further confirmed by time-averaged results. The flame behavior changes significantly when the flapping forcing is turned on. Figure 3 shows the results for the cases AF 60 and AF 75 in which the flapping forcing causes that the toroidal vortices are tilted with respect to the axial direction. Note that this effect is noticeably stronger for AF 75 . For AF 60 the tilting of the vortical rings is overwhelmed by its very strong amplification by the axial forcing term, as will be discussed in the next sections. In both the cases, however, the rings have tendency to move alternately to the left and right side of the domain and in effect the flames with the flapping excitation become wider in the 'x'-direction and narrower in the 'z'-direction. In Tyliszczak A. 19 it was shown that for low Reynolds number (Re = 4000) the flame can even bifurcate (i.e. split into two separate branches), however here this phenomenon does not occur. For the present case with Re = 23600 and relatively high level of the inlet turbulence intensity the toroidal structures are destroyed before reaching a bifurcation point that usually is located at y/D ≈ 5 33 . Figure 4a shows axial velocity spectrum computed based on a velocity signal recored in the axis at y/D = 3.0 for the case without the excitation. It can be seen that there are no distinct peaks in the spectrum that could be related to the preferred mode frequency or the parring process. Apparently, as can be seen in Figure 2a, the level of turbulence imposed at the inlet prevents formation of strong, well defined vortical structures, which periodic occurrence would certainly manifest also in the spectrum. Instead, only a rise of the fluctuations amplitudes at a broadband range of the frequencies is observed, which is centered around St = 0.6. The excitation at St a = 0.60 was chosen to match exactly this value, whereas the excitation at St a = 0.30 corresponds to its sub-harmonic at which the parring process could exist. Figure 4b shows the velocity spectra for the cases with the axial forcing only. Distinct peaks corresponding to the excitation frequencies are readily seen as they are definitely larger than the background level. The cases A 45 and A 75 do not show anything exceptional. The peaks related to their basic frequencies are virtually the only ones visible. Further downstream they become wider and lower (not presented) and it seems as there are no additional phenomena created by the forcing at these frequencies. The results for A 30 and A 60 are significantly different. In the former case (green line in Figure 4b) the excitation causes intensified velocity fluctuations not only at the basic frequency but also at its harmonic St = 0.60. The most striking difference is, however, for the case A 60 for which the whole bench of highly energetic harmonics is found. They appear as the results of interactions between subsequent vortices. These interactions take place through the wavy shape elongated rib structures (see Figure 2c) that connects the vortex rings. One could expect that existence of harmonics causes intensified mixing at small scales that speeds up the ignition process. The velocity spectra obtained for the cases with the axial and flapping forcing acting together show similar features, thought, the harmonics are much weaker. The spectra for the cases with the only flapping excitation turned on are not significantly different from the 'no forcing' case as in the axis locations close to the inlet the impact of the flapping should not be pronounced by definition.

Impact of the excitation on the lift-off height
The lift-off height (L h ) of the flames is estimated based on the time-averaged results. The time-averaging procedure started when the flames were fully developed and it lasted for the time period at least 300D/U, which was found sufficient to obtain well convergent statistics. Figure 5 shows the contours of the time averaged OH mass fraction and temperature in the central  It can be seen that both the shapes of the flames and L h are dependent on the type of the excitation. The L h was measured as the lowest point in the domain where the temperature or Y OH exceeded a given threshold. In all the cases these locations occur not in the flame axis but a few diameters off-axis. The L h predicted based on the temperature criterion is slightly smaller than using the OH mass fraction criterion, however, the differences are not very significant (∆L h < 1.0D). Worth noting is that both threshold lines predict very similar behaviour. For the case A 60 in the central part of the flame these lines are almost straight and their inclinations to the flame axis depend on the forcing frequencies (not shown). When the flapping forcing is turned on the threshold lines become rounded. Figure 6 shows dependence of the L h on the forcing frequencies for all analyzed cases. For the case without the excitation L h = 9.7D (L h = 10D in the experiment 20 ), which visibly differ from L h found applying the excitation. In these cases, depending on the forcing frequency and excitation type L h changes in the range 7.9D − 11.2D and is the lowest when the combination of axial and flapping forcing is applied. Turning on only the flapping excitation mode causes that L h reaches the maximum at St a = 0.45 after which it continuously decreases as the effect of intensified mixing of smaller spatio-temporal flow scales caused by higher forcing frequencies. For the axial excitation, both acting solely or in the combination with the flapping mode, L h behaves in a different manner. It reaches the maximum at St a = 0.75 but first it suddenly drops down for the forcing frequencies St a = 0.45 and St a = 0.60 for the axial-flapping and axial modes, respectively. The occurrence of these minima is related to the appearance of the high frequency harmonics in the spectra in Figure 4a. They have similar impact on the mixing as the increase of the forcing frequency of the flapping excitation.

Impact of the excitation on time-averaged results
It could be observed in Figure 5 that the excitation changes not only the flames positions but also visibly influences on their shapes and spreading angles. Compared to the cases with the axial excitation the flapping forcing makes the flame significantly wider. Figures 7 and 8 shows the profiles of time-averaged temperature along 'x' and 'z' directions compared with the experimental data at y/D = 1, 3, 8, 11. First, it should be noted that the present results agree well with the measurements. The location where the fuel ignites in the shear layer (y/D ≈ 11) and the near axis temperature distributions are correctly captured. Worth noting is the fact that in the experiment the co-flow temperature was biased by the 3% error 20 that could lead  to ±15 K difference compared to the assumed T cf = 1045 K. As observed by Navarro-Martinez and Kronenburg 35 the co-flow temperature has substantial impact on the flame stabilization height. Nevertheless, it seems that ±15 K error in T cf has not significant impact on the present results. Regarding the excited flames it can be readily seen that close to the inlet all profiles (also for the case without excitation) are very similar and divergences start to be seen only downstream. The temperatures are definitively larger in the axis for all excited cases, whereas for the cases AF their maxima move radially towards larger x/D locations. In general, the axial excitation causes faster temperature rise, while the flapping forcing makes the temperature more uniform along the 'x'-direction. The profiles along the 'z'-direction presented in Figure 8 shows that further from the nozzle the flapping excitation can lead to 30% narrowing of the flame in respect to the 'x'-direction.

Conclusions
The paper presented the LES studies of the H 2 /N 2 flame excited with the axial and flapping forcing. Correctness of the results was confirmed by comparison with available experimental data for the unexcited case. It was found that the lift-off height of the flame, its size and shape can be altered in a wide range depending on the type of the excitation and its frequency. Compared to the unexcited case the lift-off height can be increased or decreased. Its minimum value has been found for the case with the combination of both axial and flapping forcing at the frequency close to the center of the broadband frequency range regarded as the preferred one in the unexcited configuration. The impact of the flapping forcing manifested through a widening of the flame in the flapping direction. It was shown that the excitation can be used in advanced research on optimal control strategies.