Synchronization of thermoacoustic modes in sequential combustors

Sequential combustion constitutes a major technological step-change for gas turbines applications. This design provides higher operational ﬂexibility, lower emissions and higher efﬁciency compared to today’s conventional architectures. Like any constant pressure combustion system, sequential combustors can undergo thermoacoustic instabilities. These instabilities potentially lead to high-amplitude acoustic limit cycles, which shorten the engine components’ lifetime and therefore reduce their reliability and availability. In case of a sequential system, the two ﬂames are mutually coupled via acoustic and entropy waves. This additional inter-stages interaction markedly complicates the already challenging problem of thermoacoustic instabilities. As a result, new and unexplored system dynamics are possible. In this work, experimental data from our generic sequential combustor are presented. The system exhibits many different distinctive dynamics, as function of the operation parameters and of the combustor arrangement. This paper investigates a particular bifurcation, where two thermoa-coustic modes synchronize their self-sustained oscillations over a range of operating conditions. A low-order model of this thermoacoustic bifurcation is proposed. This consists of two coupled stochastically driven non-linear oscillators, and is able to reproduce the peculiar dynamics associated with this synchronization phenomenon. The model aids in understanding what the physical mechanisms that play a key role in the unsteady combustor physics are. In particular, it highlights the role of entropy waves, which are a signiﬁcant driver of thermoacoustic instabilities in this sequential setup. This research helps to lay the foundations for understanding the thermoacoustic instabilities in sequential combustion systems.


INTRODUCTION
In recent years, the global energy market has undergone significant changes.Both governments and customers are asking for efficient and non-polluting energy converters, and these requests are having an impact on the global energy mix.
It is well known that the increase in energy sourced from renewables raises up many technical issues, like energy storage or grid management [1][2][3].At the same time, traditional sources still play a major role.Amongst these, the combustion of fossil fuels covers 80% of the world's primary energy demands, and it is expected to maintain this role for the upcoming decades [4,5].This does not mean that the way combustion is performed will not evolve.On the contrary, the need of reducing greenhouse gases and pollutant emissions and increasing the efficiency and profitability of combustion systems drives combustion research.Gas-Turbine-based platforms are receiving increasing interest [6].First, gas turbines are naturally suited to being integrated in a network with high concentrations of renewables: thanks to their fast start-up, loading and de-loading capabilities, gas turbines can compensate the intermittency of wind and solar plants.In addition, they emit essentially no particulates, very low levels of NO x , and less than half of the CO 2 that coal-and oil-based power and propulsion systems emit.Moreover, gas turbines are fuel flexible: they can run with conventional gaseous and liquid fuels as well as with bio-and syn(thetic)-gas on demand.Over the last five years, new combustor architectures have been developed for large gas turbines and most of the main players on the market patented various axially-staged combustion concepts [7][8][9][10].All these architectures share the concept of sequential combustion: the fuel is burnt in two stages, with two flames operating at the same pressure.This allows for an increased operational range with high efficiency, ultra-low emissions and high fuel flexibility.The generic sequential combustor used in this work shares some features of the GT36 combustor [11] and is presented in the following section.

ETH Sequential Combustor
A side view of the ETH sequential combustor with some key dimensions is shown in fig. 1.This sequential combustor is operated in this study with an overall thermal power ranging from 50 to 70 kW.The air is injected through an electric heating module (upper left in the figure).The fuel for the first stage can be injected either directly into the  plenum, to achieve a full pre-mixing with air, or, like in the present study, from the fuel lance.In this case, which is closer to industrial applications, the fuel is injected from small holes in the lance downstream from the axial swirler and it is technically premixed with air inside the annular gap of the burner.Combustion takes place in the first combustion chamber, composed of a 100 mm long quartz cylinder with inner-diameter 70 mm and square cross section module of 62×62×250 mm.The tip of the lance acts as a bluffbody for flame anchoring.The vitiated gases produced by the first stage are diluted in the square section module with fresh air and then advected in the sequential burner.Here the diluted gases are mixed with the second stage fuel, injected downstream of a vortex generator.The second stage flame stabilizes in the sequential chamber, which has a square section of 88×88mm and quartz walls to provide optical access.The exhaust gases flow into an open chimney.The test rig is equipped with several measurement devices.For the present study, a water cooled Kistler piezoelectric pressure sensor type 211B2 was used to record the acoustic pressure in the first stage chamber.Additionally, two synchronized Photron high-speed camera were used to acquire images of the two flames.A sequence of the flame motion is shown on the right of fig. 1.This test rig can be used to study the fundamental physics of sequential combustors.This work focuses on thermoacoustic instability.This is a well-known problem that affects constant-pressure combustion systems and is far from being solved, even in traditional setups.Thermoacoustic instabilities originate due to a feedback interaction between acoustic pressure and the unsteady heat release rate from the flame.Under particular conditions that are difficult to predict in real engines [12], strong acoustic limit cycles can establish in the combustor.Since it is not possible to operate at these conditions to preserve the mechanical integrity of the machine, thermoacoustic instabilities limit the operational range of the turbine, or they force components to be re-designed and/or damping devices to be adopted.In the case of sequential combustion, the complexity of this problem is further increased.In this architecture, the two flames are mutually coupled, and they both contribute to the overall system stability.The second flame is particularly sensitive to temperature fluctuations, which can be generated from the first stage when the combustion is unsteady.This is especially the case in the technically premixed configuration, where coherent equivalence ratio oscillations can be part of the thermoacoustic feedback.It is expected that these entropy waves play a major role in sequential combustion thermoacoustics.In this study, experimental results collected from this sequential test rig are presented.A thermoacoustic instability involving two coupled modes is scrutinized.The dynamics of multiple modes in thermoacoustics has been studied in different contexts.For instance, in [13] two modes are alternatively active in the combustor, with stochastic transitions from one to the other.In [14] the flame response to acoustic forcing in the presence of self-sustained oscillations at another frequency is evaluated.In [15], the possibility of coexistence and coupling of two closely-spaced modes is studied.In [16], the effect of the phase of the heat release on one pair of degenerate azimuthal acoustic modes having the same frequency is analysed.The effect of synchronization phenomena between heat release and thermoacoustic field was assessed in [17].In the present case, two coexisting longitudinal modes are present.The acoustic signature in the chamber changes according to the second stage power, and therefore the system undergoes a bifurcation.The presence of an entropy-waves related mechanism as a driver of the instability is inferred.A low-order model of the two thermoacoustic modes, which is able to mimic and explain the observed dynamics is proposed.

THERMOACOUSTIC BIFURCATION
Combustion systems, like any dynamical system, are subject to bifurcations.A bifurcation is a qualitative change of the system output, when a control parameter is varied.In the present study, the effect of the second stage power on the system stability is investigated.In this study, the power of the first stage was set to Q1 =35kW, with an equivalence ratio of 0.86.The main air mass flow ( ṁA = 14g/s) was pre-heated to 70 • C, and the fuel ( ṁF1 = 0.7g/s) injected in technically premixed mode.The dilution air mass flow was ṁDA = 10g/s.The second stage power Q2 was varied, in steps of 1kW, from 15 to 35kW ( ṁF2 = [0.3;0.7]g/s), to change the total thermal power delivered.
The pressure signals were recorded for two minutes, waiting long enough to reach a stationary condition every time a new second stage thermal power was fixed.The results of this mapping are presented in fig. 2. All the operating points feature non-deterministic dynamics: the random excitation exerted by turbulence creates a probability distribution for the thermoacoustic state of the system.It is therefore interesting to examine the probability distribution of the acoustic pressure.Panel (a) shows a contour plot of P(p; Q2 ), the probability density of p as a function of the second stage power Q2 .It can be observed that the acoustic level does not change monotonically with Q2 .At 15 kW, the acoustic level is low, then it increases for intermediate powers, and finally decreases for high second stage powers.The evolution of the signal spectra when Q2 is changed was then assessed.The power spectral density of the signal recorded at three example points ( Q2 = 15; 27; 35 kW) is shown in panel (b) of fig. 2. Two thermoacoustic modes are simultaneously active in the chamber.The modes shape can be calculated via the acoustic solver AVSP, and the results are plotted in panel (d).This solver exploits a cell-vertex finite volume method of order 2. The temperature field was determined via several thermocouples placed at different axial locations on the test rig.All the inlets are designed to provide a stiff injection of air/fuel, so they are modeled as fully reflective boundaries.The outlet is equipped with an orifice plate and then discharges in a chimney.This condition was reproduced adding a hemispherical large volume in the computational domain after the outlet (not displayed in the figure) with non-reflecting boundaries.Mode I has an eigenfrequency of f I =70 Hz and exhibits high amplitudes in the plenum and in the first stage mixer, which roughly corresponds to a quarter wave resonance of the combustor.Mode II oscillates at f II =136 Hz and exhibit a shape that may be associated with a three quarter wave resonance.The first harmonic of Mode I's eigenfrequency is close to the mode II's eigenfrequency, and this fact promotes coupling phenomena between the two modes.The slight eigenfrequencies shifts are due to the change of the operating condition of temperature field and thermoacoustic feedback.To investigate the modal dynamics, the signals were band-pass filtered in two bands, highlighted in panel (b), around the two eigenfrequencies.The resulting filtered signals will be referred to as p I and p II for the first and second band, respectively.The probability density functions (PDFs) and the time traces of these two filtered signals are plotted, at the three example points, in panels (a) and (c).The two modes evolve differently when the second stage power changes.In the low power range, the spectrum is dominated by mode II, which features substantial acoustic amplitude even though no limit cycle is established.In fact, a limit cycle is characterised by a PDF with two maxima (see for instance P(p I ), shown in the second inset of panel (a)).Observing the p II time trace at the point Q2 =15kW, it emerges that it features sudden bursts of the acoustic level, followed by "quiet" periods.This behaviour is typical of a system subject to multiplicative noise [18], but this aspect will not be investigated in this work.When Q2 increases, up to 20kW, mode I slowly rises and the overall acoustic level increases.Up to this point, since mode I is linearly stable and therefore does not generate harmonics, the two modes are decoupled.At 20 kW, the overall acoustic level starts to decrease slightly and then increases again up to the maximum, which is reached at Q2 =27kW.Looking at the spectra in this range, it is possible to infer that the two modes can couple for Q2 = [20; 27] kW.In this range, since the first mode gets unstable, a new peak at the first harmonic 2 f I appear, aside from the one at f II .The increase of the overall acoustic level is due to the destabilization of mode I, which feeds with its first harmonic the second mode.As an example, the first mode exhibits a developed limit cycle at 27 kW (see again P(p I ) in the second inset of panel (a), which now features two maxima corresponding to the limit cycle amplitude).At this point, the PDF P(p II ) is unimodal but spreads to very high amplitudes.This is because the filter band includes both mode II and the harmonic of mode I.For powers higher than Q2 =27kW, the two modes decrease in intensity at the same time and they stabilize after Q2 =33kW.It is clear from this preliminary analysis that the interaction between the the observed dynamics is the result of the non trivial interaction between the two thermoacoustic modes.Further insights will be provided via the analysis carried out in the next section.

TWO-FLAMES INTERACTION
To get a further insight into the particular thermoacoustic dynamics, the filtered signals p I and p II are evaluated with a joint probability density function.This expresses the probability density of having a certain pair of values (p (i) II ), and therefore it shows how the two harmonic components of the acoustic signals are related.The results of this analysis, for six operating conditions (including the three that have been already presented), are given in the bottom row of fig. 3. It is possible to observe the evolution of the two modes instability, which is in line with what was discussed in the previous section.At low power the system is linearly stable, but the mode II features significant excursions from the rest position.Then mode I starts to appear, and a twomode limit cycle is present at 27kW: the system oscillates regularly around zero.For high power, the joint PDF shrinks towards the origin: the modes are stable again.The peculiar aspect highlighted by this plot is the oscillation "pattern" which is extremely regular.Starting from the case 19 kW, up to the most clear example of the "heart-shaped" PDF of the 27kW case, all the PDFs show asymmetries.This means that the system visits regularly and with a certain progression particular states, and so that a synchronisation of the two modes is present.To verify this, a further analysis on the oscillation phase was carried out.The single-mode thermoacoustic oscillations are slowly-modulated in amplitude and phase [19], and they can be described by: In the upper row of fig. 3 the joint PDF between the oscillation phases φ I and φ II is presented.In the lowest power case, this joint PDF is quite spread out and uniform, even though some preferential states are present.A much more narrow PDF emerges while increasing the power, together with the destabilization of mode I.In the 27kW case, three well defined high-probability density bands appear.With the help of a portion of the time trace at this operating condition, we can then explain the creation of the "heart-shaped" joint PDF for the modal acoustic pressures.The signals p I and p II are evolving together, always retracing a four-step cycle (the corresponding flames shapes were shown in fig.1): A) p I is at its local minimum (φ I = ±π), while p II is at its local maximum (φ II = 0); B) When p I crosses zero (φ I = −π/2), p II is around its minimum (φ II = ±π); C) A maximum of p I always corresponds to a maximum for p II (φ I = φ II = 0); D) p I crosses zero again (φ I = +π/2) while p II is around its minimum (φ II = ±π).Therefore the systems always follows a similar trajectory: the distribution around it is created, again, by the turbulence-induced dynamical noise, which acts both on the modes amplitudes and phases.This synchronization is present in all the operating conditions, including the very stable case of Q2 =35kW (see the joint PDF P(φ I , φ II )).However, a regular pattern is created when the oscillations significantly overcome the background noise level.The possibility of having such a synchronization is due to the fact that the second mode has a frequency that is approximately the double of the first mode's frequency (there are two cycles of p II for each period of p I ).But the question of what locks the phase difference to these specific value is still open.A possible explanation can be found in the analysis of the heat release rate signals of the two flames q 1 and q 2 .These signals can be obtained by integrating the luminosity of the high speed cameras frames over the field of view, and they approximate1 the flame heat release rate fluctuations q i (t) = V i q(x,t)dx, where V i is the volume occupied by the flame i, and q(x,t) is the overall heat release rate distribution.The corresponding spectra at the three example points are plotted in fig.4a, together with the corresponding acoustic pressure spectrum.The spectra of the heat release rate are very similar for the two flames.The two spectra also match with the pressure one.This was expected since the heat release rate oscillations are a source term for the acoustic pressure, as prescribed by the Helmholtz equation (subsonic case): where p and u are the acoustic pressure and velocity fluctuations transforms, s the Laplace variable, x the position, c the local speed of sound, γ the specific heat ratio, q the Laplace transform of the heat release rate fluctuation distribution q(x,t), n the outward normal to the boundary and Z the acoustic impedance.At each operating condition, the cross correlation χ(τ) = +∞ −∞ q 1 (t)q 2 (t + τ)dt between the heat release rate signals was evaluated.The results are plotted in fig.4b.Since the signals have a dominant component of frequency f m at each operating condition, the cross-correlation results in a modulated periodic function of the time lag τ, with period 1/ f m .In all the cases, the signals present maximum correlation for a non-zero lag.This means that the heat release rate oscillations of one flame are delayed compared to the corresponding ones in the other flames.In detail, the maximum χ is obtained for the three example cases Q2 = 15; 27; 35 kW for τ =-37.3;+11; +49.6 ms, respectively.In all of the cases, this delay is of the same order of magnitude of the convective time between the two flames.This indicates that this coupling mechanism is probably driven by entropy waves, a peculiar feature of sequential combustion system.The presence of other convective disturbances, like coherent vortical structures in the flow, was excluded by analyzing the results of a LES simulation performed in a parallel work.This time delay can set the phase difference between the two modes, as it will be shown in the next section using a surrogate model of the system.

LOW-ORDER MODEL
Low-order modelling of thermoacoustic instabilities can be done at different levels, either for predicting them, using for example thermoacoustic network models [21][22][23], or to capture the physics with simplified models that can be used for parameter identification [19,[24][25][26].In this work, a low-order phenomenological model of the thermoacoustic coupling was developed to investigate the effect of time delay between the two stages.The aim is to verify whether it is possible to qualitatively reproduce the bifurcation and to obtain a two-modes synchronized limit cycle.The acoustic pressure in the combustor is given by the superposition of several acoustic modes: where Ψ i (x) is the mode shape at the location x and η i (t) the evolution in time of the mode amplitude.The dynamics of the modes amplitudes is ruled by a set of coupled stochastic differential equations [24,27]: where ω = 2π f i is the angular frequency of the i th natural acoustic mode, and ξ i (t) is the additive stochastic forcing accounting for turbulence-induced processes weighted by the mode shape.The nonlinear function g i includes the effects of acoustic damping mechanisms and coherent2 heat release rate fluctuations.It generally depends on all the other modes amplitude and their first derivative: the coupling among modes is set by this function.
In the present case, the system dynamics is dominated by two acoustic modes, therefore we will include two modes only in our model.Regarding the non-linear function and the coupling mechanism the following aspects were considered.The damping mechanisms can often be assumed to be linear in the acoustic level.This contribution can be taken out of the function g i , and the mode equation rewritten as ηi + α i ηi + ω 2 i η i = g i (η j , η j ) + ξ i , where α is the linear damping coefficient.The RHS of this last equation represents all the terms related to the heat release rate oscillation and, rewriting in Laplace coordinates and considering again (2) and (3) yields: ( This equation states that the relation between the total heat release rate fluctuation and the mode amplitude is given by the transfer function h i .The pressure measured at the sensor position x s is then, in line with (3), p(x s ,t) The linearity of this dependency can be proven by analyzing the spectra of pressure and heat release rate, plotted for the point Q2 =27kW in fig.5a.Since p is generated by the source q (see ( 2)), if their relation were not linear, peaks having frequencies that are multiples of the frequencies at which peaks occur in q would have appeared in the p spectrum.
To complete the model, it is necessary to specify the heat release terms ξ i and g i .Regarding the first, it is possible to model this stochastic forcing as a white noise source [28].
The second term represents the flame response to acoustic perturbations.This is a non-linear and frequency-dependent function.In reality the two flames can respond to many other different perturbations, for instance to the equivalence ratio or the temperature fluctuations.In the present case we simply imagine expanding the total flame response function with a Taylor series up to third order, and we assume that it only depends on the acoustic pressure signal at the flame location (x 1 or x 2 ).In many works dealing with thermoacoustic coupling low-order modeling, the quadratic term in the Taylor expansion is absent, and this leads to a Van der Pol (VDP) oscillator equation for the mode dynamics (for instance in [29][30][31]).
In the present case, however, there are clear indications that this term should be included.Observing again the spectra in fig.2, it was already pointed out how the mode I generates, when unstable, its first harmonic at f = 2 f I .This could not be reproduced without a quadratic term in the model flame response function (FRF): the modal acoustic pressure is p I = A I (t) cos(ω I t + ϕ I ), and therefore a term p 2 I is needed to generate cos(2ω I t + 2ϕ I ) components.The presence of a quadratic term implies that the flame response is not symmetric to positive and negative acoustic pressure perturbations, as shown in fig.5b in comparison with the VDP-type FRF.Finally, to account for the presence of the two flames, we included in the model two distinct FRF b 1 and b 2 .Each FRF is made of two sub-functions, representing the response of the flame to acoustic perturbations at one of the two eigenfrequencies of the system : b f ,m (p f ,m ), m = I, II.The output of one of the two FRFs is delayed by τ D , to take into account the convective delay between the two flames.The model, including all the mentioned contributions, is presented in fig.6.It is important to highlight the role of the gain blocks, which weight each contribution by the mode shape at the location of interest.For example, the signal measured at a certain location (in the experiments, at the pressure probe location x s ) is the result of a weighting of the modes' amplitudes η m on the mode shape at the considered location Ψ x s m = Ψ m (x s ).The flame response functions take as inputs the state value at the location of the considered flame p f ,m = Ψ x f m η m .Their output is weighted again (compare to (2)) with the flame volume V f , the specific heat ratio γ and a normalization factor Λ f [19], all lumped in the terms Γ x f m = (γ − 1)/V f Λ f .What forces the acoustic blocks, finally, is the result of the weighting of the two flames' total output with the mode shape at each flame location.This is obtained by expanding (4) (compare to [19]): where q c f ,m is the coherent part of the heat release rate fluctuation of flame f at the frequency of mode m and q t f the turbulence-induced heat release rate fluctuation of flame f .To reduce the parameter space dimension, however, one can simplify this complex diagram by absorbing the gains into the functions' constants.For this purpose, one can define two new sets of functions: the already mentioned m .The resulting reduced model is depicted as a block diagram in fig.7, and the complete set of equations in time domain is: where ξ is the white noise mimicking the effect of turbulence (assumed to be the same for all the modes and flames) and p τ D m the signal p m delayed of t = τ D .As discussed, we substitute the functions B f ,m (p m ) with their expansion into a Taylor series: where β f ,m , λ f ,m and κ f ,m are the coefficients of the expansion (different for any component of the two FRF), with κ f ,m > 0 and λ f ,m either > 0 if the flame is assumed to react strongly when acoustic pressure is positive or < 0 if the response is higher when acoustic pressure is negative.The coefficient β f ,m can take any value and the phase between acoustic pressure and heat release rate can be accounted for here, for example by adopting an n-τ model [32].This set of equations was simulated in Simulink R , trying to reproduce the dynamics observed in the three experimental example points Q2 =15; 27; 35kW.The parameters of the model were tuned according to the following considerations.First, the time delay τ D was taken to be equal to the lag that gives the maximum cross-correlation χ of experimental heat release rate, and the frequency of the oscillators were chosen to be equal to the ones of the two experimental data modes.The mode shape weighting factors Ψ x m were set with the values computed with the acoustic solver AVSP (results already shown in fig.2c).The linear parameters α m , β f ,m set the linear stability of the mode, which has a linear growth rate of ν m = (α m + ∑ m β f ,m )/2.Therefore, their values are tuned to reproduce the linear stability observed in the experiments.The other constants are chosen to have oscillation amplitudes comparable with the experimental ones.The coefficients λ f ,I are necessary to link the two modes: if they are set to zero, the synchronization is not established.For simplicity, and to be able to set positive and negative delays using the same model, all the parameters are assumed to be equal over the flames, but different over the modes: a 1,m = a 2,m , a f ,I = a f ,II , ∀ parameter a.The results are presented in fig.8.The color of the boxes refers to the experimental case that the model imitates.The resulting joint PDFs of modal amplitude P(p I , p II ) and oscillation phases P(φ I , φ II ) are very similar to the ones extracted from experimental data.To assess the dependence on the initial conditions, several simulation were run, with different seeds for the random number generator of the white noise block and different initial conditions for the system of equation (7).The phase between the two oscillators locks in all of the cases to the same value, and the obtained PDFs are identical.At the same time, the desired phase synchronization is only obtained for some values of the time delay (including the ones selected to match experimental results).For other values of τ D , in fact, the phase locks to other values and the (p I , p II ) portrait changes shape.On the contrary, if the time delay is kept constant, the portrait does not change shape if any other model's parameter is varied.This fact indicates that the phase synchronization in real data might be imposed by this entropy-wave-induced time delay, which is a peculiar feature of sequential combustion system.The developed low-order model was conceived as a means to study the fundamental components of this particular setup.However, further developments are required to be able to use it as a predictive or system identification tool.In particular, additional time delays can be present in thermoacoustic systems.For instance, it is expected to have a convective delay of the swirl number fluctuations, of the equivalence ratio fluctuations, and second stage dilution air and fuel oscillations.These contributions can play a fundamental role in the definition of the dynamics of a thermoacoustic system.6 CONCLUSIONS An experimental study on thermoacoustic instabilities in a lab-scale sequential combustor was presented in this paper.The work reported how the system loses its stability for a range of values of the second stage thermal power, exhibiting a large-amplitude acoustic limit cycle.We highlighted the particular dynamics associated with this instability: two modes are active in the combustor and they interact with each other, giving rise to a synchronized bi-modal oscillation.We introduced a low-order model, composed of two coupled non-linear oscillators, which is able to reproduce qualitatively the combustor dynamics at different operative points.A key aspect of this model is its non-linear response function, which mimics the response of the flame to acoustic perturbation, and which is a cubic featuring a quadratic term.Another particular feature of the model is the presence of time delay in the feedback response.This component reproduces the delay observed between the heat release rate time traces relative to the two flame.Examining its order of magnitude, it is likely that this phenomenon is associated with the convection of entropy waves between the two stages, a very peculiar feature of sequential combustion systems.This phenomenological model was used to provide a possible explanation for the observed behavior of the combustor.Since in the model the oscillation patterns are dependent on the imposed delay between the two feedbacks, it is possible that in the real system the entropy waves traveling between the stages set a certain convective delay which defines the dynamics of the sequential combustion system.The analyzed acoustic signals presented signatures of multiplicative noise, which could be taken into account in a more sophisticated version of the low-order model, together with a more generic non-linear response function.A necessary next step is the development of parameter identification tools able to identify the value of the model's parameter, in order to gain better insights into the mechanisms of sequential combustion.The general message of this work is that sequential combustion can generate, depending on the specific setup, a variety of dynamics with the most diverse features.However, in all of the cases, entropy waves coupling is expected to be one of the key drivers of the combustion dynamics.Therefore, this aspect must be considered when dealing with thermoacoustic instability modelling and prediction.

Fig. 1 .
Fig. 1.Side view of the laboratory sequential combustor used in this study with key dimensions in [mm].All the main components are indicated with arrows.All acoustic measurements presented in this work are recorded by means of the piezoelectric sensor located in the first stage.An additional sensor is placed close to the outlet.Two series of four instantaneous snapshots of the two flames are shown.These four images correspond to an acoustic cycle, later displayed in fig.3, recorded at thermal power of first and second stages equal to Q1 =30kW and Q2 =27kW.

Fig. 2 .
Fig. 2. a) P(p; Q2 ), probability density function of the acoustic pressure oscillations p as a function of the second stage thermal power Q2 .In the insets, the PDFs at three example points ( Q2 = 15; 27; 35kW) of the acoustic pressure signal band-pass filtered in the bands 40-100 Hz (p I ) and 110-170 Hz (p II ).b) Power spectral density of the acoustic pressure at the three example operative points.The shaded areas represent the bands used to filter the signal and obtain p I and p II .The labels indicate the peaks associated with the eigenfrequency of the two main thermoacoustic modes of the chamber.c) p I and p II time traces.d) The two main acoustic modes shapes, obtained with the Helmholtz solver AVSP.

Fig. 3 .
Fig. 3. Modal dynamics in the thermoacustic bifurcation.Bottom row: joint PDF P(p I , p II ) of the modal components of the acoustic signal, at six different powers of the second stage.This joint PDFs highlights the creation of regular oscillation patterns, which is further explained via the joint PDF of the oscillation phases P(φ I , φ II ), presented for the usual three operating conditions.For the Q2 =27kW case, the acoustic filtered time traces are also shown, together with the original signal p.The labels A-D indicate the four steps of the synchronised two-modes oscillation, and they are reported for reference also on the two corresponding joint PDFs.The images of the flames at these four steps were shown in fig. 1.

Fig. 4 .
Fig. 4. a) Heat release rate spectra of the two flames, compared with the acoustic one, at the three operating conditions.b) crosscorrelation χ(τ) between the two heat release rate signals of the two flames, filtered in a band (indicated in the left panels) around the dominant frequency in the specific operating condition.The maximum χ is highlighted with a cross.

Fig. 5 .
Fig.5.a) Spectra of the total heat release rate q and acoustic pressure p.All and only the peaks of the source q are featured in p.This justifies the adoption of a transfer function to model the chamber acoustics.b) Comparison between two possible cubic models for the flame describing function.While the VDP oscillator has been already adopted and validated in the context of thermoacoustic loworder modeling, a cubic function featuring the quadratic term is necessary to reproduce the observed system dynamics.

Fig. 6 .
Fig. 6.Block diagram representing the proposed model.In the middle, the acoustic transfer functions, outputting the two modal amplitudes η m .The flames responses b f ,m are on top and on the bottom, and they take as input the acoustic modal pressure at the flame location, p f ,m = Ψ x f m η m , with Ψ x f m = Ψ m (x f ) being the mode m shape at the flame position x f .The outputs are weighted with the factor Γ x f m = (γ − 1)/V f Λ f .The flame 2 FRFs are followed by a delay block, reproducing the time lag between the two heat release rate sources.

Fig. 7 .
Fig. 7. Block diagram representing the low-order model adopted for the simulations.This model results absorbing of the gains of the model shown in fig.6 into the functions in the blocks.On the top, the two transfer functions H m (one for each eigenfrequency), reproducing the conversion of heat release rate fluctuations into acoustic pressure oscillation.On the bottom, the two Flame Response Functions, each with the two frequency components B f ,m , representing the coherent and non-linear response of the two flame to acoustic perturbations.The flame 2 FRF is followed by a delay block, reproducing the time lag between the two heat release rate sources.The output of the FRF is weighted with the mode shapes at the flames locations Ψ x f m before feeding back to the acoustic block.

Fig. 8 .
Fig. 8. Results of the Simulink R simulations of the low-order model (7), for three different sets of system parameters.In the panels, the joint PDFs for the two filtered output signals and phases.Cases (a), (b), and (c) reproduces qualitatively the dynamics observed at the three operating points Q2 = 15, 27, 35kW, respectively.