Nonlinear Dynamical Features of Vortex-Acoustic Lock-On in a Backward-Facing Step Combustor

Joel Varun Vasanth 1

Indian Institute of Technology, Madras Chennai 600036, India

Satyanarayanan Chakravarthy Indian Institute of Technology, Madras Chennai 600036, India

ABSTRACT Unstable shear layers characteristic of step combustors display strong vortex shedding and play a dominant role in combustion instability. Past non-premixed flame experiments with Reynolds number (Re) as the bifurcation parameter showed that during instability, the frequency of acoustic oscillations locks-on to the natural vortex shedding mode for a range of Re. In the present work, we study the dynamical features of this type of lock-on. We employ a vortex model to bring out the nonlinearity in the flame response to vortical perturbations. This is coupled with a Galerkin model for the acoustic field. Regimes of aperiodic oscillations, low-amplitude lock-on and high-amplitude lock-on are predicted in agreement with experimental trends. A transition to instability during lock-on via a supercritical Hopf bifurcation is observed with an increase in limit cycle amplitudes. Synchronisation between the acoustic pressure and vortex-driven heat release perturbation is studied using recurrence analysis. The system transitions from asynchronous to phase synchronisation in the low-amplitude lock-on regime, characterised by an increasing degree of phase correlation. In the high-amplitude lock-on regime, a state of generalised synchronisation exists, where in addition to the phase, the amplitudes also show strong correlation.

1. INTRODUCTION

Flow fields past flame-holding devices in industrial and aircraft burners consist of spatially developing shear layers with regions of reverse flow that are inherently hydrodynamically unstable [ 1 ]. The resulting vortex shedding is a major factor contributing to combustion instability [ 2 ].

A unique feature of vortex-dominated combustion instability is vortex-acoustic lock-on. This has been experimentally observed by varying parameters of the flow such as the air flow Reynolds number, Re , in several flow configurations - a non-premixed backward-facing step (BFS) combustor [ 3 ], a partially premixed blu ff -body combustor [ 4 ] and recently a blu ff -body combustor

1 joelvarunvasanth@gmail.com

a slaty. inter.noise 21-24 AUGUST SCOTTISH EVENT CAMPUS O ¥, ? GLASGOW

with syngas fuel mixtures [ 5 ]. In [ 3 ], it was observed that, at a certain Re in the non-premixed BFS combustor, the dominant mode switches from f a to f v and remains locked on to f v with further increase in Re, accompanied by a nonlinear rise in p ′ RMS amplitudes. Similar observations of lock-on were made in [ 4 ] by varying combustor length, both air and fuel flow rates, and injection and blu ff -body locations. Baraiya et al [ 5 ] observed lock-on with the shear layer mode and the third and sixth harmonic of the longitudinal duct acoustic modes.

Modelling approaches must include a vortex-based flame response in addition to an acoustic flame response to incorporate the competing e ff ects of either mode during lock-on. Vortex kicked-oscillator models [ 6 ], [ 7 ], [ 8 ]. They assume an unsteady heat release rate proportional to the circulation of a growing vortex at the location of the step edge. However, a clear-cut lock-on regime has not been predicted since a framework for dynamical evolution of the flame to capture flame-vortex interactions is lacking.

Many studies on lock-on in the context of externally forced thermoacoustic systems have recently been performed [ 9 ], [ 10 ], [ 11 ]. These studies also characterise the nonlinear dynamics of lock-on. General results show a lock-on as a synchronisation of the self-excited thermoacoustic oscillation with the forcing frequency as the forcing amplitude is increased, via a quasiperiodic route. Low order models to mimic such dynamics were proposed, such as the forced van der Pol oscillator [ 11 ] and the Landau-Stuart model [ 12 ]. However, the modes considered in the above studies were all acoustic in nature and contributions from hydrodynamically self-excited oscillations were not dealt with.

The forcing of hydrodynamically self-excited modes were studied by [ 13 ] and [ 14 ] in unconfined jet di ff usion flames, where lock-in of the hydrodynamic mode to the forcing frequency was observed. These studies, including those above, involve a unidirectional coupling from the applied forcing to the acoustic or hydrodynamic modes. What is required from the point of view of lock-on, is a bidirectional coupling. Further, these studies included scenarios where the vortex shedding oscillation frequency locks-on to the acoustics and not the other way around. Thus, developing a model to capture this latter type of lock-on and study the dynamical characteristics is of interest to the current work.

We develop individual submodels for each of the three participating processes – flame, vortex shedding and acoustics, to enable their simultaneous evolution in a fully three-way coupled manner. For the non-premixed flame, we extend the use of a two-dimensional thermo-di ff usive model in a BFS combustor to include the e ff ect of the vorticity field, relaxing the assumption of uniform flow. For the 2D flow, we perform a Schwartz-Christo ff el conformal mapping of the step geometry and represent the velocities using complex potentials. This overcomes the localised vortex burning assumption in kicked-oscillator models and introduces the Arrhenius-type and convective nonlinearities of the acoustic driving as imposed on the base flow. Along with the Galerkin method for the acoustics, the resulting coupled simulations are compared with experimental results. Synchronisation between the acoustic pressure and heat release perturbation is then studied using phase space reconstruction.

2. MODELLING

2.1. Flow Modelling We assume an two-dimensional, inviscid, laminar flow in a backward facing step (BFS) combustor, whose dimensions are adopted from [ 3 ] and illustrated in Figure 1 (a). A uniform flow of velocity ˜ U

enters at the inlet. A Schwartz-Christo ff el transformation (SCT) is used to map this geometry from the physical plane ˜ z = ˜ x + j ˜ y onto the upper half of a transformed plane ˜ ξ = ˜ ξ x + ˜ ξ y . In general, the transformation is defined as

d ˜ z d ξ = ˜ K Y

α i

π − 1 , (1)

i ( ξ − p i )

where the RHS multiples over the ξ vertices of the geometry p i that subtend corresponding angles α i in the physical plane ˜ z and ˜ K is a scaling constant which may be complex. Here, dimensional quantities are represented by a tilde. Complex potentials ˜ w are described in the ξ plane, and their complex velocities ˜ σ in the ˜ z plane are then obtained from the relation

˜ σ = d ˜ w

d ˜ z (2)

Following this, the physical velocity field an be extracted as ˜ u = ℜ ( ˜ σ ) , ˜ v = ℑ ( ˜ σ ). The SCT as applied to the BFS is displayed in Figure 1 and is given by

1 − η − 1 √ a log √ a + η √ a − η

s

d ˜ z d ξ = ˜ K π

" log 1 + η

a − ξ 1 − ξ (3)

# , η =

Figure 1: (a) Two-dimensional BFS geometry (black lines with red labels for dimensions) as adapted from [ 3 ] in the physical ˜ z plane showing fuel and oxidiser inlets and direction of flow (blue labels) (b) Schwartz-Christo ff el transformation of the geometry in the ξ plane showing the point source (blue arrows with filled arrowheads) at the ξ origin B’C’.

Next, we decompose the flow field ˜ u into three components – a base flow ˜u B devoid of vortices, the velocity field induced by vortices ˜u V and the acoustic velocity ˜ u ′ . The last component is obtained as explained in section 3 . To obtain ˜u B , we define a potential source of strength ˜ m = ˜ U ( ˜ k + α ˜ h ) /π at the ξ origin given as ˜ w B = ˜ m log ξ .

We now find ˜ K and a by substituting 3 in 1 and evaluating it at ξ = 0 and ∞ where σ = ˜ U and ˜ U (˜ k + α ˜ h ) / ˜ H respectively. Doing this, we get ˜ K = ˜ H /π and a = ˜ H 2 / (˜ k + α ˜ h ) 2 . To model ˜u V , we adopt the Rankine vortex model with a vortex strength ˜ κ at a location ξ V . In order to ensure no through-flow boundary condition across the walls ξ y = 0, it is necessary to introduce an image vortex in the ξ y < 0 region with equal strength but opposite sign. The potential is then

˜ w V = − j ˜ κ log ( ξ − ξ V ) + j ˜ κ log  ξ − ¯ ξ V  (4)

with the overbar denoting complex conjugate. The vortex circulation ˜ Γ relates to ˜ κ as ˜ κ = ˜ Γ / 2 π . In the presence of N V vortices and their images, the combined complex velocity field ˜ σ = ˜ σ B + ˜ σ V in the ˜ z plane is

(b) éplane by A nN Bic! & oO E' F

N V X

κ i ξ − ξ V , i κ i ξ − ¯ ξ V , i (5)

˜ σ = ˜ m

ξ + j

i

A rotational core of radius ˜ r A centered at each vortex centre ξ = ξ V , i is adopted in order to remove the irrotational singularities in the Rankine vortex model. The vorticity ˜ Ω is thus constant at all locations within the core and is equal to 2˜ κ/ ˜ r 2 A . The velocities at grid points within ˜ r ∈ [0 , ˜ r A ] are replaced by tangential velocities ˜ Ω ˜ r . This constitutes the velocity fields within the vortex core.

Next, we proceed to model the dynamics of vortex growth and shedding. The growth of circulation in a vortex that is formed behind the step as a result of shear layer separation is found as the integral of vorticity generated per unit time over the thickness of the boundary layer at the step edge. The contribution of the remaining boundaries to the creation of vorticity and growth of circulation is assumed to be absent. The accumulation of circulation leading to the roll-up of a vortex is modelled using a quasi-steady hypothesis from [ 6 ] in terms of the velocity at the outer edge of the boundary layer ˜ u s

d ˜ Γ

d ˜ t = ˜ u s ( t ) 2

2 (6)

where ˜ u s ( t ) = ˜ U + ˜ u ′ ( t ), ˜ u ′ being the acoustic velocity at the step. When the strength of the growing vortex reaches a critical value ˜ Γ cr = ˜ u s ˜ h / 2 S t , it sheds. The Strouhal number S t is natural frequency ˜ f v of vortex shedding scaled as S t = ˜ f v ˜ h / ˜ U . We set the St to correspond to the passage frequency of the large-coherent structures formed after merging of smaller vortices at the step as a result of collective interaction at instability. Thus, we assume the absence of vortex-vortex interactions such as merging of previously shed vortices. The St is taken as the value defining the vortex shedding frequencies pertaining to the BFS geometry in [ 3 ] which is the geometry adopted in the present model. These St are experimentally obtained, similar to [ 6 ] and [ 7 ]. The shed vortices are advected at the free stream velocity ˜u B , and the total velocity field computed by Eq. 5 is evaluated on a rectangular 2D grid, and used as the velocity for the flame model described next.

2.2. Flame Model We consider a non-premixed combustion configuration as shown in Figure 1 , with the fuel inlet being a fraction α of the step height. We use the thermo-di ff usive model [ 15 ], as derived from the Navier-Stokes equations for reacting flow, with the main assumption being that the density is a constant. This can be used to model premixed and di ff usion flames and flame with varying premixedness [ 16 ]. We modify this to include the two-dimensionality of the BFS flow described in the previous section, thus relaxing the uniform flow assumption considered in previous adopted models. We adopt finite rate chemistry and solve for the mass fraction of fuel ˜ Y , oxidiser ˜ X and temperature ˜ T in the reacting field. The convective nonlinearity overcomes the linear relationship in previous kicked oscillator models [ 6 ], [ 7 ]. The density of the mixture, the heat capacities of both species, kinematic viscosity of the mixture, and the molecular di ff usivities are assumed constant and the Lewis number is taken as unity for all species.

Non-dimensionalisation is performed with reference quantities ˜ T re f = ˜ Q / ˜ C p , where ˜ Q is the fuel heating value and ˜ C p is the mixture specific heat; ˜ Y re f = ˜ α F = ν F ˜ W F , ˜ X re f = α O = ν O ˜ W O , where ν and ˜ W are the stoichiometric coe ffi cient and molecular weight of each reactant respectively. The pressure is scaled by ˜ ρ 0 ˜ c 2 0 /γ . We non-dimensionalise the spatial coordinate with the step height as x c = ˜ x / ˜ h and time with the flow time scale as t = ˜ t / ( ˜ h / ˜ U ). We thus generate a set of three 2D equations for the solution variable Φ = [ X Y T ] T . We then rewrite the dimensionless equations as

D Φ

Dt = 1

Pe ∇ 2 Φ + ⃗ s DaXY exp ( − θ/ T ) (7)

We adopt flux boundary conditions for Φ at the inlets and no flux boundary conditions at the walls and step. The far upstream values for fuel and oxidiser mass fractions used in the flux boundary conditions - X i and Y i are provided in Table 1 . The integrated heat release rate per unit length is given as where ⃗ s = [1 1 − 1]. Using the mixture di ff usivity D, mixture density ˜ ρ 0 , activation energy ˜ E a , universal gas constant R u and pre-exponential factor B, we form the Peclet number Pe = ˜ U ˜ h / ˜ D , Damkohler number Da = ˜ B α F α O ˜ h / ˜ ρ 0 ˜ U and scaled activation energy θ = ˜ C p ˜ E a / ˜ R u ˜ Q . The velocity in the convective term on the LHS of (2.14) is the total field u .

q x = Z 1

− 1 DaXYexp  − θ

 dy c (8)

T

and the total heat release rate :

q = Z 1

0 q x dx c (9)

In order to calculate the integrated heat release rate perturbation per unit length q ′ x , we use a base state heat release q x , B defined as that obtained from a ‘base state flame’ in the base flow u B excluding vortices. Then, q ′ x = q x − q x , B [ 15 ]. The total heat release rate perturbation q ′ is then calculated as q ′ = R 1

0 q ′ x dx c .

2.3. Acoustic Model We assume one-dimensional longitudinal mode acoustics due to the non-observance of transverse modes in the given geometry for the range of flow conditions considered. We assume that e ff ects on the acoustics due to non-uniformity in the temperature field and finite Mach number e ff ects due to the mean and oscillating flows are absent. The dimensional linearized acoustic momentum and energy equations are then written as

∂ ˜ p ′

∂ ˜ u ′

∂ ˜ x = 0 , ∂ ˜ p ′

∂ ˜ t + ˜ ρ 0 ˜ c 2 0 ∂ ˜ u ′

∂ ˜ t + 1

∂ ˜ x = ( γ − 1) ˜ q ′ (10)

˜ ρ 0

where ˜ c 0 is the speed of sound and γ is the ratio of specific heats. The non-dimensionalisation is performed as x a = ˜ x / ˜ L a . Here we take ˜ L a ≡ ˜ L , the length of the combustor duct as shown in Figure 1 . Time is scaled with the flow time scale ˜ h / ˜ U used thus far in the flow and flame models above, since the time marching is performed along with the flame equations at the flow time step. The dimensionless equations then contain a time scale factor τ = ˜ L c / M ˜ L a , M being the mean inlet flow Mach number. A perfectly open-open duct is assumed for the model. Then, using the Galerkin method, u ′ and p ′ are written in terms of the corresponding natural mode shapes

N M X

N M X

j d η j

j = 1 cos ( j π x a ) η j , p ′ = − γ M

1

u ′ =

dt sin ( j π x a ) . (11)

πτ

j = 1

where N M is the number of Galerkin modes η j . An ad hoc mode-wise damping ζ j is adopted [ 6 ]. If ω j is the angular frequency of the j th mode,

! 0 . 5   (12)

 c 1 ω j ω 1 + c 2

ω 1

ζ j = 1

2 π

ω j

where we take the damping coe ffi cients c 1 = 0 . 27 , c 2 = 0 . 03 as measured in a BFS [ 17 ]. Since the flow is incompressible, the Mach number is low and we assume an acoustically compact combustion

zone with the location of acoustic forcing at the step x f = ˜ x f / ˜ L a . It is also assumed that any time delay incurred in pressure waves travelling from the location of vortex burning to the step is negligible. Thus, for the j th mode [ 15 ]:

d 2 η j

dt 2 + 2 ζ j k j τ d η j

dt +  τ k j  2 η j = − 2 τ k j

M q ′ sin  j π x f  , k j = j π (13)

which is then solved with the initial conditions

η 1 = 0 . 1 , η j = 0 . 0 , j = 2 , 3 , 4 , . . . N M ; d η j

dt = 0 , j = 1 , 2 , 3 , . . . N M (14)

The equations for the combustion and acoustic fields are integrated in time by a 4th order Runge- Kutta (RK4) scheme. The convective components of the flame equations 7 are discretized with a first order upwind scheme and the di ff usive components are discretized with a second order central di ff erence scheme. Grid independence tests were performed and a grid size of 121x241 is adopted, along with a time step of 10 − 4 . A list of remaining parameters used in the model are given in table 1 .

Table 1: List of parameters used in the current study along with their values.

Parameter Value Parameter Value Parameter Value Parameter Value

ν F 1 θ 1.07 Q 45 MJ / kg ˜ k 0.03m

ν O 9.52 ˜ Y 1 Pe 10-30 a 0.13

W F 16 g / mol ˜ ρ F 0.656kg / m 3 Da 2.9-9.9 × 10 1 0 β 2

W O 28.9g / mol ˜ E a 36.5kcal ˜ h 0.03m N M 10

Re 10000-60000 ϕ 0.043-0.26 ˜ Y i 1 T i 0.05

3. RESULTS

We first present the results of the vortex growth and shedding model in non-reacting flow and without the influence of acoustic coupling in section 3.1. In section 3.2, we proceed to use the present model to predict vortex-acoustic lock-on over a range of Re in the BFS geometry in [ 3 ]. In the following section 3.3, the synchronisation analysis is performed on the data from the model. ,

3.1. Flow Field Results These results are in Figure 2 . First, we show the irrotational streamlines in part of the geometry close to the step in Figure 2 (a) as a result of the base flow u B devoid of vortices. Large velocities are present at the step edge O due to the irrotational field. Figures 2(b)-(d) show the growing vortex embedded in the flow field u B + u V governed by Eqs. (2.9)-(2.11). A region of reverse flow is formed, starting from the step and ending at the reattachment point. Due to the growing circulation, the size of the region of reverse flow becomes larger, as can be seen by the advancement of the reattachment point in figures 2 (b)-(c). At t = 6.67 = 1 / St, the shedding condition is reached, the growing vortex is shed and is advected further away from the step at the base flow velocity u B .

Figure 2: (a) Irrotational streamlines pertaining to velocity field without vortices u B in the dimensionless z plane. (b) Vortex growth at t = 5, shedding at (c) t = 8, growth of a new vortex after shedding of a previous fully grown vortex at (d) t = 11 . 5.

3.2. Vortex-acoustic lock-on We now proceed to solve the present model over a range of Re as in [ 3 ] as described in table 1 . The acoustic pressures p ′ are plotted in figure 3 over frames (a)-(h) for various Re tried. The lowest Re (frame (a)), shows the absence of any e ff ective self-excitation. The p’ fluctuations remnant from the initial condition decay rapidly due to the inherent damping in the model. At a slightly higher Re (frame (b)), shows some sustained p’ oscillations at a non-negligible amplitude. However the amplitudes are not large enough to indicate a strong feedback loop. In the experimental case, this low Re range 10000-20000 is akin to stable operation with low amplitude broadband turbulent combustion noise. Hence, we refer to this as the ‘stable regime’.

With a further increase in Re (frames (c)-(f)) a steep rise in peak-to-peak amplitudes is observed. At Re = 35000-45000, limit cycle oscillations are observed at a well-defined frequency. Limit cycle behaviour persists over this extended range of Re. This corresponds to an ‘unstable regime’. These high amplitude tonal self-sustained oscillations in the unstable regime are characteristic of combustion instability. Further increasing the Re, a drop in amplitudes is evident (frames (g)-(h)).

The spectral contents of the above p ′ evolution are now studied. The dominant frequencies over the entire range of Re tried, are plotted in Figure 4 (a), and their corresponding amplitudes are plotted in Figure 4 (b). Plotted for comparison with the present results are the corresponding experimental and LES results from [ 18 ]. The model predicts low amplitude p’ oscillations at the natural duct acoustic frequency ˜ f a 105 Hz from in the range Re = 10k to 25k, corresponding to the modes at stable operation in the experiments and LES. At these Re , a large disparity in the time scales of natural acoustic and vortex shedding oscillations prevents an e ff ective feedback loop, responsible for the observed low p’ amplitudes. At Re = 27500, a jump in the dominant mode of p’ oscillation from ˜ f a to the frequency of forcing, i.e. the fundamental vortex shedding mode ˜ f v is observed, corresponding to the onset of instability. The Re of onset compares fairly well with the experimental and LES values of 30000 and 26000 respectively, considering the simplicity of the present model. The dominant mode of p ′ oscillation continues to occur at ˜ f v up to an Re of 50000. The accompanying nonlinear rise in amplitudes with Re is evidence that this regime is a vortex-acoustic lock-on.

The trends in the amplitude variation (Figure 4 (b)) with increase in Re are captured well by the model. A sustained low amplitude trend prior to lock-on followed by a rise in amplitudes at the

2 wi (a) (b) Tl 2 oxidisex SNH 0 fel) oe (>—--§$_ =, 0 0 1 Xe 2 3 0 1 X, 2 3 (©) 2 @ Fl th —_—e—E (|) BReESeeS 0 1 2 3 0 1x 2 3

Figure 3: Evolution of acoustic pressure p ′ as a result of three-way coupling between flame, flow and acoustic fields for the range of Re as obtained from the current model.

onset of lock-on Re = 30000 is observed. From the amplitude plot, it is clear that the nature of the bifurcation to instability as the parameter Re is increased is a supercritical Hopf bifurcation, given the smooth increase in amplitude at each Re without any abrupt jumps. The amplitudes peak at Re = 40000 when ˜ f v = ˜ f a . This is true for the experimental and LES data as well. However, the amplitude comparison must be treated as largely qualitative. Even the LES results [ 18 ] overestimate the amplitudes because of inadequate modelling of the acoustic damping. In the present model, this is further compounded by the simplicity of the flow model adopted wherein the turbulent fluctuations and the energy loss from the large scale flow to turbulence is not included. Further work to include turbulent e ff ects is currently being pursued.

3.3. Synchronisation In order to participate in syncrhonisation, the oscillators must be self-sustained, i.e., they must have a unique energy source, damping, and characteristic frequency [ 19 ]. In the thermoacoustic system, the acoustics and the flow modulated heat release constitute a pair of such oscillators, and so we use the signals p’ and q’ and extract their phases. Methods of analysis based on recurrence plots (RP) are common for oscillators that are non-phase coherent. The RP method as defined in [ 20 ] is elucidated below.

Recurrence is a fundamental characteristic of a dynamical system. Recurrence plots are binary matrices R i , j that represent whether a state x j at time instant j recurs to an ϵ -neighbourhood of its value at a former time i. It is represented as R i , j = Θ ( ϵ || x i − x j || ), where Θ is the Heaviside function. Based on this definition, we use a metric called probability of recurrence or generalized auto-correlation function P ( τ ) to quantify synchronisation. P ( τ ) is defined as the probability of such a recurrence after a certain time delay τ .

0.01 F (a) Re = 10000 r (b) Re = 20000 0 jw WW evernnrnnnnnnnnnennne » -0.01 | L 5 1 1 m 1 1 G 0.01 | (c) Re = 30000 L (a) Re = 35000 oy 2 0 NW g 3 0.01 T _, , | _ 001 (¢) Re = 40000 (f) Re= 45000 & 2 9 = 0.01 = : ; ; ; : Z 0.01 + (g) Re = 50000 L (h) Re = 60000 : -0.01 + r 0 30 60 90 120 0 30 60 «690 )=6120 150 Nondimensional Time

Figure 4: (a) Dominant mode frequencies and (b) corresponding amplitudes of the present model plotted along with the experimental [ 3 ] and LES [ 18 ] data at similar flow conditions. For reference, inclined dotted lines in (a) are plotted corresponding to modes of vortex shedding and horizontal dotted lines for the natural acoustic modes, as labelled.

The phase space is reconstructed via the time-delay embedding method where a time delay T D and embedding dimension E d are to be determined à priori. The T D is obtained using the concept of average mutual information, where the optimum T D is that which corresponds to the first minimum of mutual information. If N r = N − ( E d − 1) T D is the total number of reconstructed vectors for N data points, and ϵ is a recurrence threshold limiting the nearness of the states,

N r − τ X

P ( τ ) = 1 N r − τ

i = 1 R i , i + τ (15)

The above method enables us to identify phase and generalised synchronisation (PS and GS), usually detected in combustion systems [ 8 ]. PS occurs when the phases and frequencies of the oscillators are locked but their amplitudes remain uncorrelated. A locking of the positions of local maxima of P ( τ ) for q ′ and p ′ (or P q ′ ( τ ) and P p ′ ( τ ) respectively) indicates PS since both signals recur at the same intervals. A much stronger synchronisation is when the amplitudes are correlated in addition to the phases and frequencies. When the oscillators synchronise in such a way, their trajectories have a functional relationship and they are said to be in GS. This features a matching of the magnitudes in addition to the positions of P q ′ ( τ ) and ( P p ′ ( τ )).

We further elucidate on what synchronisation properties are possessed during extended lock-on period presented in this work. The P ( τ ) at Re’s corresponding to stable condition, onset of lock-on, and unstable condition of limit cycle oscillations at maximum amplitude are extracted in figure 5 . The P q ′ ( τ ) in all three cases show well-pronounced local maxima. This is due to the modulation of the heat release by the vortex shedding at a definitive frequency. The peaks occur regularly at the distinct frequencies of their respective signals. For Re = 10000 (Figure 5 (a)), we see the number of peaks in the p ′ plot is larger than that of the q ′ since the natural acoustic mode at this Re is much larger. Yet, the system is apparently in PS, since some P q ′ ( τ ) peaks coincide with P p ′ ( τ ) at regular intervals. However, it is actually not in PS, as we shall show by the definition of an additional parameter in the next paragraph. At the onset of lock-on

Experiment » Computation ° Present model 5 250 (eq) sopmydury Jeuonemndwuos ‘eournedxg o 2 S 7150 10 0 100 200 300 400 500 600x100 escoeoeceoe S$ssSs8 88585 Sua FRA (ea) sopmydury jepoyy wasorg 2 = rt * : s 2 af oP Ma 41s a: © : a : EI oc 6 1s sy 18 ry 18 + r=) Ss a ry Ss a J J S esescese $5 SS & a A a 4 (zy) Aouonbary yueurm0g Re Re

Figure 5: Probability of recurrence showing the nature of synchronisation between q ′ and p ′ at Re (a) in stable regime, (b) at onset of lock-on and (c) at limit cycle behaviour. The recurrence threshold ϵ is taken as 10% of the maximum phase space diameter, T D = 500ms and E D = 20.

at Re = 30000 (figure 5 (b)), we see a perfect correspondence in the positions of the peaks, implying PS. The driving is not su ffi ciently strong to cause a GS, as can be seen from the mismatch in their peak magnitudes. The transition to GS occurs at Re = 50000 (figure 5 (c)) where we observe the largest dominant mode amplitudes (figure 4 (b)). This observation can be attributed to the e ff ective coupling of the two oscillators due to the low phase di ff erences and the closeness of their frequencies, enabling a one-to-one functional relationship between their trajectories, indicating strong driving.

12 1 0.8 5 0.6 0.4 0.2 0 0 0.04 0.08 0.12 0.16 0.2 T (s) (a) Re = 10000 0 0.02 0.04 0.06 T (8) (b) Re = 30000 0 0.01 0.02 0.03 0.04 7 (s) (c) Re = 50000

As mentioned above, we turn our attention to case Re = 10000 again to clarify if the trajectories are in PS or not. For this purpose, we calculate the coe ffi cient of probability of recurrence (CPR) throughout the Re range. The CPR is defined as CPR = ⟨ ¯ P p ′ ( τ ) ¯ P q ′ ( τ ) ⟩ / ( χ p ′ χ q ′ ), where ⟨·⟩ denotes the mean, overbar denotes the mean subtracted values, and χ is the standard deviation of P ( τ ) [ 20 ]. It serves as a quantitative indicator of PS by calculating the linear correlation between the P ( τ ) data. A pair of trajectories showing a high CPR 1 tend to recur at similar times scales, implying they are in PS.

Figure 6: Coe ffi cient of Probability of Recurrence (CPR) showing the degree of phase synchronisation across the range of Re.

CPR os 06 O48 02 10000 20000 30000 40000 50000 60000 Re

We plot the CPR in figure 6 . A low CPR value of 0.3 at Re = 10000 implies p ′ and q ′ are indeed asynchronous, although their frequencies are equal. This is simply because of their phases not being fully locked due to insu ffi cient driving amplitudes, leading to the stable regime. However, this condition ceases at the onset of lock-on, where CPR shows abrupt jumps to higher values, subsequent to which the high CPR values prevail through the lock-on range, exhibiting that they are phase locked. The extended Re range of high CPRs elucidates how synchronisation is fundamentally di ff erent from resonance. The q ′ − p ′ interaction being two-way coupled, shows a locking on when the forcing frequency is over a wide range around the natural frequency. Resonance is a one-way coupled phenomenon where high oscillation amplitudes are observed only when natural frequency detuning is small. The two-way coupling as modelled here enables synchronization by a simultaneous evolution of the two processes. Even though at the start of the simulations the detuning is considerable, as the flow acoustics and flame evolve simultaneously, the acoustically generated large vortex roll-up makes the q ′ strong enough to cause a phase locking. This results in the high correlation in P ( τ ) trajectories.

4. CONCLUSION

The present work attempts to capture the onset and frequency range of lock-on using a reduced- order model. We introduce a three-way coupled model, modelling the vortex shedding, reacting field, and acoustic field, based on their individual physics. The SCT-based potential flow formulation for growth and shedding of vortices is used as the flow model. The thermal-di ff usive model is used to model the flame and the acoustic equation is modelled with the Galerkin method.

The acoustic pressure evolution at each Re showed a variation in the amplitudes from low to high amplitude limit cycle oscillations. The acoustics oscillates at its own natural mode frequency during stable operation. At a certain Re, a jump in this dominant mode frequency to the vortex shedding natural mode is observed. The Re of onset of lock in is well predicted and so is the range of frequencies.

Syncrhonisation between the p ′ and q ′ signals is analysed by recurrence analysis. The signals transition from being asynchronous prior to lock-on, to phase synchronous and finally generally synchronised after the onset of lockon. Once phase synchronisation is achieved, the simultaneous evolution draws the phases closer, facilitating a stronger forcing and thus a steep amplitude rise. This is the fundamental mechanism of the lock-on phenomenon.

REFERENCES

[1] P Huerre. Local and Global Instabilities in Spatially Developing Flows. Annual Review of Fluid Mechanics , 22:470–537, 1990.

[2] K.C. Schadow and E. Gutmark. Combustion Instability Related to Vortex Shedding in Dump Combustors and their Passive Control. Progress in Energy and Combustion Science , 18:117–132, 1992.

[3] S. R. Chakravarthy, R. Sivakumar, and O. J. Shreenivasan. Vortex-acoustic lock-on in blu ff -body and backward-facing step combustors. Sadhana - Academy Proceedings in Engineering Sciences , 32(1- 2):145–154, 2007.

[4] R. Sivakumar and S. R. Chakravarthy. Experimental Investigation of the Acoustic Field in a Blu ff -Body Combustor. International Journal of Aeroacoustics , 7(3-4):267–299, 2008.

[5] Nikhil A. Baraiya and S. R. Chakravarthy. Excitation of high frequency thermoacoustic oscillations by syngas in a non-premixed blu ff body combustor. International Journal of Hydrogen Energy , 44(29):15598–15609, 2019.

[6] Konstantin I Matveev and F E C Culick. A model for combustion instability involving vortex shedding. Combustion Science and Technology , 175, 2003.

[7] Vineeth Nair and R. I. Sujith. A reduced-order model for the onset of combustion instability: Physical mechanisms for intermittency and precursors. Proceedings of the Combustion Institute , 35(3):3193–3200, 2015.

[8] Samadhan A. Pawar, Akshay Seshadri, Vishnu R. Unni, and R. I. Sujith. Thermoacoustic instability as mutual synchronization between the acoustic field of the confinement and turbulent reactive flow. Journal of Fluid Mechanics , 827:664–693, 2017.

[9] Saravanan Balusamy, Larry K.B. Li, Zhiyi Han, Matthew P. Juniper, and Simone Hochgreb. Nonlinear dynamics of a self-excited thermoacoustic system subjected to acoustic forcing. Proceedings of the Combustion Institute , 35(3):3229–3236, 2015.

[10] Yu Guan, Meenatchidevi Murugesan, and Larry K.B. Li. Strange nonchaotic and chaotic attractors in a self-excited thermoacoustic oscillator subjected to external periodic forcing. Chaos , 28(9), 2018.

[11] Karthik Kashinath, Larry K.B. Li, and Matthew P. Juniper. Forced synchronization of periodic and aperiodic thermoacoustic oscillations: Lock-in, bifurcations and open-loop control. Journal of Fluid Mechanics , 838:690–714, 2018.

[12] Hiroaki Hyodo and Tetsushi Biwa. Phase-locking and suppression states observed in forced synchronization of thermoacoustic oscillator. Journal of the Physical Society of Japan , 87(3):1–5, 2018.

[13] Matthew P. Juniper, Larry K.B. Li, and Joseph W. Nichols. Forcing of self-excited round jet di ff usion flames. Proceedings of the Combustion Institute , 32 I(1):1191–1198, 2009.

[14] Larry K.B. Li and Matthew P. Juniper. Phase trapping and slipping in a forced hydrodynamically self- excited jet. Journal of Fluid Mechanics , 735:R5, 2013.

[15] M. Tyagi, S. R. Chakravarthy, and R. I. Sujith. Unsteady combustion response of a ducted non-premixed flame and acoustic coupling. Combustion Theory and Modelling , 11(2):205–226, 2007.

[16] David S Bhatt and S R Chakravarthy. Nonlinear dynamical behaviour of intrinsic thermal-di ff usive oscillations of laminar flames with varying premixedness. Combustion and Flame , 159(6):2115–2125, 2012.

[17] James D. Sterling. Nonlinear Analysis and Modelling of Combustion Instabilities in a Laboratory Combustor. Combustion Science and Technology , 89(1-4):167–179, 1993.

[18] Ashwin Kannan, Balaji Chellappan, and Satyanarayanan Chakravarthy. Flame-acoustic coupling of combustion instability in a non-premixed backward-facing step combustor: the role of acoustic-Reynolds stress. Combustion Theory and Modelling , 11(4):658–682, 2016.

[19] Arkady Pikovsky, Michael Rosenblum, and Juergen Kurths. Synchronization: A Universal Concept in Nonlinear Sciences . Cambridge, Potsdam, 1 edition, 2001.

[20] Norbert Marwan, M. Carmen Romano, Marco Thiel, and Jürgen Kurths. Recurrence plots for the analysis of complex systems. Physics Reports , 438(5-6):237–329, 2007.