# Numerical investigation of azimuthal thermoacoustic instability in a gas
turbine model combustor

Zhi X. Chen<EMAIL_ADDRESS>Nedunchezhian Swaminathan  Marek Mazur
Nicholas A. Worth  Guangyu Zhang  Lei Li  State Key Laboratory of Turbulence
and Complex Systems, Aeronautics and Astronautics, College of Engineering,
Peking University, Beijing 100871, P.R. China  AI for Science Institute
(AISI), Beijing, 100080, P.R. China  Department of Engineering, University of
Cambridge, Trumpington Street, Cambridge CB2 1PZ, UK  CORIA-UMR 6614 Normandie
Université, CNRS-Université et INSA de Rouen, Campus Universitaire du
Madrillet, Saint Etienne du Rouvray, France  Department of Energy and Process
Engineering, Norwegian University of Science and Technology, N-7491 Trondheim,
Norway  Research Institute of Aero-Engine, Beihang University, Beijing 100191,
P.R. China

###### Abstract

Self-excited spinning mode azimuthal instability in an annular combustor with
non-swirling flow is investigated using large eddy simulation (LES).
Compressible Navier-Stoke equations are solved with a flamelet combustion
model to describe the subgrid chemistry$-$turbulence interactions. Two
flamelet models, with and without heat loss effects, are compared to elucidate
the non-adiabatic wall effects on the thermoacoustic instability. The
azimuthal modes are captured well by both models with only marginal
differences in the computed frequencies and amplitudes. By comparing with the
experimental measurements, the frequencies given by the LES are approximately
10% higher and the amplitudes are well predicted. Further analysis of the
experimental and LES data shows a similar dominant anti-clockwise spinning
mode, under which a good agreement is observed for the phase-averaged heat
release rate fluctuations. Dynamic mode decomposition (DMD) is applied to shed
more light on this spinning mode. The LES and experimental DMD modes
reconstructed for their azimuthal mode frequencies agree very well for the
heat release fluctuations. The DMD mode structure for the acoustic pressure
from the LES shows a considerable non-zero profile at the combustor outlet,
which could be essential for azimuthal modes to establish in this annular
combustor. Finally, a low-order modelling study was conducted using an
acoustic network combined with the flame transfer function extracted from LES.
The results show that the dominant mode is associated with the plenum showing
a first longitudinal and azimuthal mixed mode structure. By tuning the plenum
length to match the effective volume, the predicted frequency becomes very
close to the measured value.

###### keywords:

Large Eddy Simulation; Annular Combustor; Swirl flame; Thermoacoustic
oscillation; Self-excited instability.

††journal:

##  1 Introduction

Azimuthal thermoacoustic instabilities remain an unresolved issue for gas
turbine (GT) industries [1, 2]. This instability arises when the unsteady
combustion of multiple circumferentially aligned flames couples with the
acoustic resonant modes of the annular chamber, which often dominates over
other modes exhibiting higher oscillation frequencies and amplitudes [3]. As
pressure waves propagate not only in the longitudinal direction but also in
the transverse direction along the circumference, azimuthal modes involve more
complex dynamics such as flame/flame interaction [4], modal dynamics [5],
azimuthal/longitudinal mode coupling (e.g., the slanted mode) [6], etc. These
complexities make azimuthal modes very difficult to understand, predict and
control.

To simulate and predict thermoacoustic instabilities, a range of approaches
have been developed over the past several decades, from hand-calculated
original analytical [7, 8] models, then desktop running low-order network
models, to large eddy simulation (LES) now performed on supercomputers,
accompanied with a vast increase of computational cost [2]. Due to the complex
burner geometry and boundary conditions, analytical and low-order network
models still play essential roles in studying azimuthal modes, see Refs. [9,
10] for extensive reviews on these models. To reduce the computational burden,
hybrid approaches combining incompressible LES and acoustic solvers have also
been developed [11, 12], where the connection between solvers, i.e., the flame
response to acoustic perturbations is considered through a flame transfer or
describing function (FTF or FDF). The advantage of this approach is that the
flame calculation, e.g. reactive flow LES, is decoupled from the system
acoustics calculations so that it can still be performed using the
computationally much cheaper low-Mach formulations. The major difficulty,
however, is how to apply appropriate forcing to the flame when the nozzle
inlet involves complex geometry with multiple swirlers, or the perturbations
need to be applied in transverse directions for azimuthal instabilities. With
the rapid advances in high-performance parallel computing over recent years,
compressible LES for full annular lab-scale burners has now become
computationally feasible [13, 14], offering a promising alternative to
investigate azimuthal instabilities and the complex underlying physical
mechanisms for modal dynamics and structures.

Due to the high cost of running annular test-rigs, few laboratory annular
combustor experiments are available (e.g., [5, 15]). It has been shown in
these studies that the swirl-generated azimuthal mean flow can introduce
complex time- and operating condition-dependent switching among spinning,
standing and mixed modes. To better understand the modal behaviour of spinning
mode without azimuthal bulk flow, Mazur et. al [16] removed the swirlers from
the 12-burner annular rig studied in [5, 4] and showed well-characterised
spinning modes occurring for high amplitude instabilities over a wide range of
operating conditions. Therefore, these experiments provide an ideal test case
for LES of azimuthal instabilities. The main objective of this study is to
simulate the non-swirling annular premixed combustor of Mazur et. al [16] and
explore the LES capabilities of capturing spinning modes observed in the
experiments. Moreover, previous LES studies have been performed over
relatively short run times (tens of ms, much shorter than experiments). The
present study also investigates the effect of running for extended times in
order to evaluate the time varying mode development.

Amongst the available subgrid combustion models, flamelet approaches are
preferable for complex geometries because of their robustness and low
computational cost. However, compared to models that compute chemistry
explicitly (not necessarily refer to computing chemistry using resolved
quantities known as no model or quasi-laminar chemistry, but also in thickened
flame and partially stirred reactor models, where local temperature and
pressure effects are inherently considered in reaction rate closure),
tabulation models require additional non-trivial modelling treatment to
account for heat loss effects, which were shown to have a strong impact on the
FTF/FDF [12]. More recent studies by Kraus et al. [17] and Agostinelli et al.
[18] employed advanced conjugated heat transfer method to account for the two-
way coupling between the flame LES and burner walls. They showed that the heat
loss effects are important for the swirl burners considered and need to be
taken into account in the LES in order to correctly capture the thermoacoustic
instability. Thus, a non-adiabatic flamelet approach is adopted from Massey et
al. [19] for this work and the heat loss effects on the instability simulation
are investigated.

This paper is organised as follows. The target experimental configuration is
described in Section 2, followed by the LES modelling framework and numerical
setup in Section 3. The LES results are discussed in Section 4 and the low-
order modelling is presented in Section 5. The conclusions are summarised in
Section 6.

##  2 Experimental test case

The target annular combustor was first experimentally studied at Cambridge by
Worth and Dawson [5, 4] to characterise and understand azimuthal modes in
laboratory-scale swirled burners at atmospheric pressure. By varying the
number of burners (i.e., flame-to-flame spacing) and operating flow
conditions, various azimuthal modal dynamics were observed involving complex
coupling and interplay among different physical processes, most of which still
remain to be understood. To simplify the system and thus better understand the
individual processes, this combustor was recently modified by Mazur et. al
[16] at NTNU by removing the swirlers in the 12 flames configuration (largest
flame spacing). As such, azimuthal modes could be established in the absence
of azimuthal mean flow and flame$-$flame interactions, in a geometrically
simpler flow.

Fig. 1: Mid-plane cut view of the experimental configuration [16] (left) and
typical LES snapshot undergoing ACW spinning mode (indicated by white arrows).
Zoomed view of the mesh near bluff-body is also shown.

The experimental geometry has been described in detail in the previous work
[5], and is only summarised briefly here together with important changes. A
schematic of the modified burner is shown in Fig. 1 with all the dimensions in
mm. A perfectly premixed ethylene-air mixture is passed through a cylindrical
plenum chamber, before being divided around a central bluff body into 12
injector tubes. Each injector tube holds a central rod, to which a conical
bluff is mounted with a diameter of 13 mm. In the present work, a bulk flow
velocity of 24.5 m/s and a equivalence ratio of 0.95 are employed. Further to
Ref. [16], the geometry was slightly modified in order to more closely
replicate the simulation boundary conditions. In particular, all flow
conditioning elements were removed from the plenum and a high level of
uniformity was demonstrated suggesting that the conditioning provided by the
plenum is sufficient. An additional plate was mounted inside the inner annular
wall cylinder (marked as ‘blocking plate’ in Fig. 1) in order to prevent the
entrainment of ambient air through the central opening during the experiment.
Measurements [16] were performed using dynamic pressure sensors ($f_{s}=51.2$
kHz, sampled over $t_{s}=10$ s) at 3 equi-spaced angular locations, as shown
in Fig. 1. Pressure ports were flush-mounted in the injector tubes, at a
distance $45$ mm upstream of the dump plane. Heat release rate imaging is
performed using a high speed intensified camera ($f_{c}=10$ kHz, sampled over
$t_{c}=2$ s, spatial resolution 3.2 pixels/mm).

The computational domain includes the full geometry of the combustor including
the upstream plenum, 12 premixers (feeding tubes), combustion chamber and the
downstream atmospheric buffering domain. Note that the conical volume at the
bottom of the plenum was omitted and only the cylindrical part (marked by
green in Fig. 1) is considered to simplify the inflow conditions. However,
this simplification has an considerable impact on the predicted mode
frequencies, which will be discussed in detailed in Section 5. A typical
snapshot of LES undergoing an anti-clockwise (ACW) spinning mode is also shown
in Fig. 1. Chamber pressure contours are plotted with the arrows indicating
the direction of the travelling azimuthal wave. The oscillating flame surfaces
are marked using volume rendering for the reaction rate of progress variable.
The axial velocity contours in the bluff-body mid-plane are shown for two
opposite burners. From one flame to another along the annulus, the flame
length varies in response to the acoustic perturbations resulting in azimuthal
heat release fluctuations, which give rise to further acoustic perturbations
in the feedback loop. To capture this phenomenon, an appropriate LES framework
considering compressibility effects is required along with correct numerical
boundary conditions for the acoustics. In addition, a robust and efficient SGS
combustion model is also essential and these are described next.

##  3 LES framework and numerical setup

The LES modelling approach used for this study involves the flamelet chemistry
of 1D freely propagating laminar premixed flames and a presumed probability
density function (PDF) to account for the subgrid scale (SGS)
turbulence$-$chemistry interactions. This approach has been previously tested
and validated extensively for a broad range of flames and configurations. In
particular, it was shown in [20, 21] that the unstable flame dynamics and a
longitudinal (Helmholtz) mode thermoacoustic instability in a dual swirl GT
burner were well captured. Thus, it is of interest here to further assess this
approach for azimuthal modes in annular combustors. The modelling framework is
briefly described below and see details in [21] and Refs. therein.

Compressible 3D Navier-Stokes equations are solved for filtered mass and
momentum in a pressure-based form for the low-Mach number flows considered in
this study:

$\displaystyle\frac{\partial\,\overline{\rho}}{\partial
t}+\nabla\cdot\left(\,\overline{\rho}\,\widetilde{\bm{U}}\,\right)$
$\displaystyle=0\,,$ (1)
$\displaystyle\frac{\partial\,\overline{\rho}\,\widetilde{\bm{U}}}{\partial
t}+\nabla\cdot\left(\,\overline{\rho}\,\widetilde{\bm{U}}\,\widetilde{\bm{U}}\,\right)$
$\displaystyle=-\nabla\,\overline{p}~{}+~{}\nabla\cdot{\bm{\tau}}_{\rm
eff}\,,$ (2)

where the modified filtered pressure $\overline{p}$ is computed by solving a
Poisson equation 111It is worth noting that here the formulation used is
different from the so-called variable-density low Mach solver, in which
pressure and density are solved in a decoupled manner and hence density
variation is only a result of varying temperature. The implementation of
$\rho$-PIMPLE algorithm in OpenFOAM-v7 is employed and it has been shown in
previous studies [21, 22] that this algorithm is able to capture aero- and
thermo-acoustic phenomena accurately.. The filtered mixture density,
$\overline{\rho}$, is obtained using the ideal-gas equation of state
considering pressure variations. The effective shear stress tensor is given by
${\bm{\tau}}_{\rm
eff}=2\overline{\rho}\left(\,\widetilde{\nu}+\nu_{t}\right)\left[\widetilde{\cal{\bm{S}}}-(\nabla\cdot\widetilde{\bm{U}})\,{\textbf{I}}/3\right]$,
where $\widetilde{\cal{\bm{S}}}$ is the filtered strain rate,
$\widetilde{\nu}$ and $\nu_{t}$ are the filtered molecular and subgrid eddy
viscosities respectively. The latter is computed using the Smagorinsky model
[23]. To describe SGS turbulent premixed combustion with heat losses, the non-
adiabatic flamelet thermochemical states are pre-tabulated using three control
variables: filtered reaction progress variable, $\widetilde{c}$, its SGS
variance, $\widetilde{c^{\prime\prime 2}}$, and total enthalpy,
$\widetilde{h}$ (sensible$+$formation), which are solved using their transport
equations:

$\displaystyle\frac{\partial\,\overline{\rho}\,\widetilde{c}}{\partial
t}+\nabla\cdot\left(\,\overline{\rho}\,\widetilde{\bm{U}}\;\widetilde{c}\,\right)$
$\displaystyle=\nabla\cdot\left(\,\overline{\rho}{\mathscr{D}}_{\rm
eff}\nabla\,\widetilde{c}\;\right)\,+\,\overline{\dot{\omega}_{c}}\;,$ (3)
$\displaystyle\frac{\partial\,\overline{\rho}\,\widetilde{c^{\prime\prime
2}}}{\partial
t}+\nabla\cdot\left(\,\overline{\rho}\,\widetilde{\bm{U}}\,\widetilde{c^{\prime\prime
2}}\,\right)=\;$
$\displaystyle\nabla\cdot\left(\,\overline{\rho}{\mathscr{D}}_{\rm
eff}\nabla\,\widetilde{c^{\prime\prime
2}}\;\right)\,-2\,\overline{\rho}\,\widetilde{\chi}_{{\rm sgs}}\,$
$\displaystyle+2\,\overline{\rho}\frac{\nu_{t}}{{\rm
Sc}_{t}}|\nabla\,\widetilde{c}\,|^{2}+2\left(\,\overline{c\,{\dot{\omega}_{c}}}-\widetilde{c}\,\overline{\dot{\omega}_{c}}\,\right)\,,$
(4) $\displaystyle\frac{\partial\,\overline{\rho}\,\widetilde{h}}{\partial
t}+\nabla\cdot\left(\,\overline{\rho}\,\widetilde{\bm{U}}\,\widetilde{h}\,\right)$
$\displaystyle=\nabla\cdot\left[\,\overline{\rho}\left(\frac{\widetilde{\nu}}{\rm
Pr}+\frac{\nu_{t}}{{\rm
Pr}_{t}}\right)\nabla\widetilde{h}\,\right]+{\frac{D\overline{p}}{Dt}}\,,$ (5)

where the laminar and turbulent Prandtl numbers, ${\rm Pr}$ and ${\rm
Pr}_{t}$, are both taken to be 0.7. The effective diffusivity
${\mathscr{D}}_{\rm eff}$ is modelled using $(\,\widetilde{D}+\nu_{t}/{\rm
Sc}_{t})$ with a turbulent Schmidt number ${\rm Sc}_{t}=0.4$ [24, 25], and
$\widetilde{D}$ is the filtered molecular diffusivity. The progress variable
is defined using the sum of CO and CO2 mass fractions and it is normalised to
vary from 0 to 1 between the fresh and burnt gases. The unclosed reaction rate
and dissipation rate terms are modelled through a revised flamelet formulation
as detailed in Refs. [21, 26]

Laminar 1D premixed flame are calculated using Cantera for a range of heat
loss levels. The full USC mechanism (C1$-$C3) [27] is used for ethylene/air
combustion at 300 K and 1 atm. The simulation results are plotted against the
experimental measurements [28] in Fig. 2.

Fig. 2: Laminar flame speed versus equivalence ratio for ethylene/air mixture
at 300 K and 1 atm. Symbols: measurements [28]; Thick line: adiabatic
simulation; Thin lines: simulations with heat losses.

Heat loss effects are implemented in a similar fashion to that used in Ref.
[19] by scaling down the heat release rate in the energy equation. As can be
seen in the figure, the computed adiabatic flame speeds agree well with the
measurements and as the heat loss increases, the flame speed drops smoothly.
Figure 3 compares the measured and computed mean flame shapes using time-
averaged OH∗ chemiluminescence from the experiments and heat release rate from
LES with non-adiabatic (NAD) and adiabatic (AD) flamelet models. Note that
since the OH∗ image is a line-of-sight visualisation, the LES results are
integrated in the transverse direction. It is seen that the flame shape and
length are well captured in the LES. In comparison with the AD model, the NAD
case shows weaker burning outer branches with extinguished flame roots near
the dump plane, which is similar to the experimental image. Bauerheim et. al
[12] argued that the lifted flame base is more prone to acoustics-induced
oscillations, which may have a impact on the flame response for azimuthal
modes. This is further discussed later in Section 4.2.

Fig. 3: Mean shape of single flame: line-of-sight integrated LES heat release
rate using adiabatic (left) and non-adiabatic (middle) flamelets, compared
with the experimental [16] OH∗ image (right).

The LES modelling framework described above is implemented using the OpenFOAM
2.3.0 libraries. The computational domain includes the upstream tubes and
plenum, and a fixed mass flow rate condition ($u^{\prime}=0$) is employed for
the inlet boundary. To mimic the open downstream boundary, a large
hemispherical domain is added and the far-field boundary is specified as a
pressure outlet ($p^{\prime}=0$). A fully unstructured grid with about 18
million ($\sim$1.5M per burner) tetrahedral cells is used and a zoom-in view
of the mesh near the bluff-body is shown in Fig. 1. Mesh sensitivity was
investigated for a single sector using a finer mesh and similar results were
found for the mean fields and also the FTF (see Appendix A). Isothermal
condition of 1000 K [14, 12] is given for the combustion chamber and bluff-
body top walls and adiabatic condition is specified for all other walls. The
Time-step is determined such that CFL $<$ 0.5 for the entire domain and the
typical size is between $2-5\times 10^{-7}$ s. A fully developed cold flow is
first simulated for about 200 ms before ignition and then there is a transient
period for azimuthal modes to settle and subsequently establish limit-cycle
oscillations. The length of this transient period ranges between tens and
hundreds of ms, highly depending on the ignition method used and this is
discussed next.

##  4 Results and discussion

###  4.1 Onset of computed azimuthal modes: effect of the ignition method

For LES of multi-sector combustors, ignition is commonly conducted by
duplicating the reactive fields of a separate single-sector simulation to
reduce the cost [14]. This differs from the experiments [16] where an ignition
source was used to light the burner from downstream and subsequently the flame
propagated upstream until stabilised around the bluff-bodies. Since the
ignition process (known as light-around) involves an azimuthal motion of the
flame, which gives rise to azimuthal pressure waves, the ignition transients
may have an impact on the start-up of azimuthal instabilities by acting as an
initial transverse forcing.

Fig. 4: Pressure time series for two different ignition cases (for adiabatic
simulations). Probes P1, P2 and P3 are marked in Fig. 1 (upper tapping
position).

To shed light on this, two separate LES cases, referred to as full ignition
(FI) and spark ignition (SI) respectively hereafter, are performed based on
the same initial fields. Figure 4 presents the computed time series of the
pressure fluctuation at probes P1, P2 and P3, which are located in the
injector tubes, 120∘ apart in the azimuthal direction and 45 mm upstream of
the dump plane. For the FI case, the entire annular chamber is ignited by
duplicating the single-sector solutions. The ignition location for SI is on
the centreline of one of the bluff-bodies and the height is the same as the
inner wall (see marked in Fig. 1). For both the FI and SI cases, a transient
evolution from ignition to fully excited azimuthal mode is observed in the
LES. Shortly after numerical ignition at $t=0.195$ s, significant longitudinal
oscillations ($>$500 Pa) appear and the three pressure signals collapse onto
one another. This period is quite long for the FI case and the large
oscillations occur several times spanning over almost 0.1 s. For the SI case,
the longitudinal oscillations reach a much higher amplitude, about twice as
that for the FI case, which consequently initiate the transition period at a
much earlier time about 0.015 s after ignition. During the transition period,
mixed longitudinal and azimuthal oscillations co-exist and the length of this
period is quite similar between the two ignition cases, although the
amplitudes are generally higher for the SI case. Finally, dominant azimuthal
modes are established in both cases showing very similar frequency and
amplitude as can be seen in the zoomed signals in Fig. 4. This ignition
sensitivity test suggests that although azimuthal instability can be
reproduced using either ignition methods, the onset of the azimuthal modes is
strongly influenced by the transient ignition processes and the spark ignition
approach followed by a natural light-around process leads to a significantly
faster establishment of the instability. This has a significant implication
for LES of azimuthal instabilities as the computational cost can be
substantially reduced.

###  4.2 Comparison with the experiments

The statistics obtained from experiments and simulations for the spinning
azimuthal modes are compared in this subsection. It should be noted that the
start-up process, i.e., from ignition until limit-cycle oscillation, was
observed to be quite long (up to almost 10 s), although the azimuthal
instability already appears in the first several hundreds of milliseconds.
During this process, the oscillation frequency remains nearly the same from
the beginning, whereas the amplitude ramps up gradually. As simulating the
entire process for 10 s is computationally impractical, the same sampling
period, 400 to 500 ms after ignition, was chosen for both the experimental and
LES pressure data to make an appropriate comparison.

Figure 5a compares the measured pressure power spectra and those obtained from
LES with the AD and NAD approaches described earlier in Section 3. It can be
seen that both models are able to capture the azimuthal instability although
the predicted frequencies are 100$-$150 Hz higher than the measured value at
around 1800 Hz. The reason for this difference is that the acoustic mode
structure dictating the excited mode lies in the upstream plenum, having one
wavelength in both the longitudinal and azimuthal directions. Therefore,
omitting the conical volume (mentioned earlier in Fig. 1) effectively resulted
in a shorter longitudinal length and hence a higher fundamental frequency.
More details on this will be discussed later in Section 5.

Fig. 5: Comparison of measured and computed (a) pressure power spectra and (b)
joint PDF of the azimuthal wave amplitudes with the spin ratio (SR) plotted in
the inset.

The amplitude in Fig. 5 is also under-predicted by both LES models by 5$-$8%.
By including the heat loss effects, the NAD model gives a slightly lower
frequency and smaller amplitude, but the overall difference is marginal
compared to the adiabatic approach. This suggests that unlike the compact
swirling flames [12], the bluff-body stabilised flames are much longer and
their response to acoustic perturbations are less affected flame root dynamics
induced by the wall heat losses. This is also consistent with previous
numerical studies on the FTF/FDF of a single bluff-body flame [29], see e.g.
Refs. [30, 31], where the non-linear flame response was well captured by
adiabatic models. As the two LES models give similar results, further analyses
hereafter are performed only for the adiabatic case for simplicity.

Fig. 6: Comparison of phase-averaged fluctuating heat release rate (FHRR)
obtained from overhead OH∗ chemiluminescence measurements (top) and line-of-
sight integration of LES results (bottom).

To compare the azimuthal mode characteristics, the measured and computed
pressure data at probes P1, P2 and P3 (see Fig. 1) are processed to obtain the
super-position of azimuthal wave amplitudes following the least-square fitting
procedure used in Refs. [16, 5]. The amplitudes of anti-clockwise and clock-
wise spinning waves are denoted as $A^{+}$ and $A^{-}$ respectively. Based on
these amplitudes, the spinning and standing modes can be well defined and the
joint PDFs of $A^{+}$ and $A^{-}$ obtained from the measurement and LES are
compared in Fig. 5b. A good agreement can be observed with both exhibiting an
evident ACW-spinning mode. While the experimental PDF shows a considerable
scatter, the LES results fall very close to the vertical axis with high PDF
values, suggesting a stronger preference in the ACW spinning direction. This
may be indicative of less noise in the LES [32, 33]. It is also to be noted
that in this symmetric burner configuration (without bulk swirl), there should
be equal probabilities for CW or ACW modes. Some level of asymmetry due to
manufacturing or other uncertainties resulted in a preferences for ACW
observed over many experimental runs. However, the geometry is perfectly
symmetric in the LES, which should prefer neither directions, and the ACW mode
observed is likely to be a coincidence for a single LES run. Using the
definition given in [15], the spin ratio can be calculated as: ${\rm
SR}=(A^{+}-A^{-})/(A^{+}+A^{-})$. The measured and computed SR are plotted in
the inset of Fig. 5b. Consistently with the joint PDFs, both the LES and
experiment give SR values close to 1 throughout the sampling period, featuring
a reasonably stable ACW spinning mode. Mode switching is not observed for this
case due to the relatively high equivalence ratio and bulk velocity, which was
shown to prefer spinning modes [16].

The strength and distribution of the fluctuating heat release rate (FHRR) and
its response to the spinning azimuthal waves are of prime importance for
understanding azimuthal instabilities. Thus, phase-averaged normalised FHRRs
for six angles ($0^{\circ}$, $60^{\circ}$, …, $300^{\circ}$) obtained from the
experimental OH∗ images and LES results are compared in Fig. 6. The maximum
heat release fluctuation in the LES is about 7%, which agrees well with the
experimental value of 6%. The previously reported [16] crescent shaped
structures around both sides of the flame is also observed here in both the
measurements and simulations. Despite the considerable noise in the
experimental data, a quite symmetric FHRR distribution is seen with half of
the burners forming a positive semicircle and the other half forming a
negative one. This structure spins as one moves along the cycle from left to
right in the figure.

The above comparisons for the pressure and heat release rate show an overall
good agreement between the experiment and simulation, and this suggests that
the spinning mode and its major characteristics are well captured in the LES.

###  4.3 DMD analysis of the azimuthal mode

Dynamic mode decomposition (DMD) [34] is a powerful technique to study
coherent spatial-temporal structures in dynamic systems. Compared to proper
orthogonal decomposition (POD), DMD generates a set of modes associated with
their own oscillation frequency, amplitude and growth/decay rates, which are
desired for a better interpretation of combustion dynamics [35]. Hence, DMD is
applied here to further investigate the spinning mode azimuthal instability.

The DMD algorithm implementation of [34] is followed to process 516 snapshots
collected from continuous measurements and simulations spanning over about 30
ms. This corresponds to approximately 50 cycles of the azimuthal mode. Note
that the experimental snapshots are 2D line-of-sight OH∗ images taken from the
top of the burner, and the LES snapshots are 3D heat release rate field for
the entire annular combustion chamber. Prior to the DMD computation, the LES
results are integrated in the streamwise direction and then projected to the
same 2D plane as the OH∗ images. The DMD mode growth rates for the peak
frequency (about 1800 and 1900 Hz for the experiment and LES respectively) are
near-zero negative values, suggesting that limit-cycle oscillations have been
established before the sampling period. It is worth remarking here that in the
DMD spectrum for the LES, two eigen-frequencies at about 900 and 1100 Hz also
showed considerable amplitudes, which are consistent with the two small peaks
seen in Fig. 5 and the low frequency content seen on the time series in Fig.
4. These are associated with the longitudinal modes of the combustor geometry.
However, both of these DMD modes have large negative growth rates suggesting
that they will eventually be damped out in a longer simulation. Nevertheless,
the azimuthal mode of interest is well established already in the present
simulation duration and hence is discussed next.

Fig. 7: Phase-shifted DMD reconstruction of the FHRR at the azimuthal mode
frequency for three phase angles. The corresponding combustion chamber
pressure mode is shown for the LES.

The DMD reconstructions of the heat release rate at the azimuthal oscillation
frequencies are shown in Fig. 7 for three phase angles. The reconstructed
pressure field (not available from the experiments) is also shown for the LES.
Note that these phases here do not necessarily correspond to those shown in
Fig. 6, because there is an arbitrary phase shift during the DMD process. It
is seen that the reconstructed modes show relatively less noise compared to
the simple phase averages, particularly for the LES. Again, good agreement is
obtained here for the DMD modes for FHRR, which further confirms that the
spinning mode behaviour is well captured in the LES. It is also seen that the
DMD modes for pressure and FHRR seem to be synchronised in the annular
chamber. Finally, the DMD mode structure for the LES domain is shown in Fig.
8. The mode structure features a first longitudinal/first azimuthal mixed mode
in the upstream plenum, a first longitudinal mode with the pressure node at
the middle of the premixer, and a quarter-wave longitudinal/first azimuthal
mixed mode in the combustion chamber. These are found to be quite similar to
that computed in [12] using a Helmholtz solver. However, evident differences
are observed here for the downstream combustion chamber outlet (top of the
outer wall), where the LES acoustic pressure retains a considerable value
compared to the $p^{\prime}=0$ condition imposed in the Helmholtz solver. This
downstream mode structure may be essential for self-excited azimuthal modes to
occur. It also implies that downstream acoustic boundary condition could have
a major role on the overall stability and azimuthal modes, and hence should be
handled appropriately in the analytical and low-order modelling approaches.

Fig. 8: DMD mode structure of the acoustic pressure reconstructed at the
azimuthal mode frequency from the LES.

##  5 Low-order modelling

In order to further understand the azimuthal mode predicted in the LES and its
dependency on the upstream plenum geometry, a low-order modelling analysis is
performed. The model combustor geometry is simplified to establish an
analytical network acoustic model, as shown in Fig. 9.

Fig. 9: Sketch of the half cross-section view of the simplified combustor. The
sketch is not to scale.

Note that the network model considers the system acoustics in the three-
dimensional space (see Ref. [36] for details) and the sketch here only shows a
half cross-sectional view of the 3D model. We consider the pressure
disturbance in the form of

$p^{\prime}=\mathrm{Re}[\hat{p}\,\mathrm{e}^{\mathrm{i}\omega t}],$ (6)

where the operator $\mathrm{Re}[\cdot]$ takes real part of the complex value.
The premixers are so thin that the waves inside can be assumed to be planar
waves, written as

$\hat{p}_{1N}(x)=P_{1N}^{+}\mathrm{e}^{\frac{-\mathrm{i}\omega
x}{\overline{c}_{1}+\overline{u}_{1N}}}+P_{1N}^{-}\mathrm{e}^{\frac{\mathrm{i}\omega(x-L_{1})}{\overline{c}_{1}-\overline{u}_{1N}}},$
(7)
$\hat{u}_{1N}(x)=\frac{1}{\overline{\rho}_{1}\overline{c}_{1}}[P_{1N}^{+}\mathrm{e}^{\frac{-\mathrm{i}\omega
x}{\overline{c}_{1}+\overline{u}_{1N}}}-P_{1N}^{-}\mathrm{e}^{\frac{\mathrm{i}\omega(x-L_{1})}{\overline{c}_{1}-\overline{u}_{1N}}}],$
(8)

where $P_{1N}^{+}$ and $P_{1N}^{-}$ denote the complex amplitudes of the
downstream and upstream propagating waves respectively, and $N$ denotes the
sequence number of the premixers. The axial coordinate $x$ is relative to the
left end of each network element (0 for plenum, 1 for premixer and 2 for
combustion chamber). The local speed of sound and mean velocity are denoted by
$\overline{c}$ and $\overline{u}$ respectively, with the subscript
corresponding to the given network element. The waves in the plenum and
combustion chamber are three-dimensional, and can be obtained from the wave
equation using common notations:

$\displaystyle\hat{p}_{j}(\mathbf{r})=\sum_{m=-\infty}^{m=+\infty}\sum_{n=1}^{+\infty}$
$\displaystyle\psi_{jm}(z_{jmn}r)\times$ (9)
$\displaystyle[P_{jmn}^{+}\mathrm{e}^{ik_{jmn}^{+}x}+P_{jmn}^{-}\mathrm{e}^{ik_{jmn}^{-}(x-L_{j})}],$
$\displaystyle\hat{u}_{j}(\mathbf{r})=-\frac{1}{\overline{\rho}_{j}}\sum_{n=1}^{+\infty}\psi_{jm}(z_{jmn}r)$
$\displaystyle[\frac{k_{jmn}^{+}P_{jmn}^{+}\mathrm{e}^{\mathrm{i}k_{jmn}^{+}x}}{\omega+k_{jmn}^{+}\overline{u}_{j}}$
(10)
$\displaystyle+\frac{k_{jmn}^{-}P_{jmn}^{-}\mathrm{e}^{\mathrm{i}k_{jmn}^{-}(x-L_{j})}}{\omega+k_{jmn}^{-}\overline{u}_{j}}],$

where $j=0,2$ and $\psi$ denotes the eigenfunction of the circular duct
($j=0$) or annular duct ($j=2$), $m$ the azimuthal mode number and $n$ the
radial mode number. The wavenumbers $k_{jmn}^{\pm}$ can be obtained by
dispersion equations [36].

Since the combustion zone is nearly isobaric, the flame can be seen as compact
and is located at the dump plane [10]. The convective terms involving
vorticity and entropy are neglected. Hence, the energy conservation across the
flame is reduced to

$\frac{\gamma}{\gamma-1}\frac{\partial}{\partial
x_{j}}(\overline{p}\hat{u}_{j})=\hat{q}.$ (11)

The unsteady heat release rate term $\hat{q}$ using the linear $n-\tau$ flame
transfer function (FTF) is obtained from the LES (full-annular adiabatic
case), with the form

$\frac{\hat{q}_{N}}{\overline{q}_{N}}=k_{f}\frac{\hat{u}_{1N}}{\overline{u}_{1N}}e^{-\mathrm{i}\omega\tau},\;x=L_{0}+L_{1},$
(12)

where the gain $k_{f}=1.1$ and the phase difference $\omega_{\tau}=1.77$. The
reference point for the velocity fluctuations was taken to be 60 mm upstream
of the dump plane to keep clear of the inlet turbulence disturbances. It is
worth remarking that the choice of this location has a quite considerable
effect on the above FTF values due the standing mode inside the premixers.
Also note that the FTF is only reliable at the self-excited frequency of the
LES results. Equation (11) is integrated over the ($r$, $\theta$) plane [37].
Then we have the closed equations that characterize the disturbances in the
model combustor

$\mathbf{X}(\omega)\mathbf{P}=\mathbf{0}.$ (13)

By solving $\det[\mathbf{X}]=0$, we have the complex eigenfrequency
$\omega=\omega_{r}+\mathrm{i}\omega_{i}$ of the system, where $\omega_{r}$
represents the natural frequency, and $-\omega_{i}$ represents the growth
rate. The positive growth rate means the system is unstable, while the
negative growth rate means stable. Refer to Ref. [36] for elaborate modelling
details.

Table 1:  Model geometrical parameters. Parameters | Symbols | Values
---|---|---
Length of plenum | $L_{0}$ | 0.2 m
Radius of plenum | $R_{0}$ | 0.106 m
Length of premixer | $L_{1}$ | 0.155 m
Cross-section area of premixer | $A_{p}$ | 0.00015 $\mathrm{m}^{2}$
Length of combustion chamber | $L_{c}$ | 0.3 m
Outer radius of combustion chamber | $R_{2}^{o}$ | 0.106 m
Inner radius of combustion chamber | $R_{2}^{i}$ | 0.0635 m
Table 2:  Frequency variation against the plenum length. $L_{0}/\mathrm{m}$ | 0.12 | 0.14 | 0.16 | 0.18 | 0.20 | 0.22 | 0.24
---|---|---|---|---|---|---|---
$f/\mathrm{Hz}$ | 3009 | 2637 | 2386 | 2115 | 1900 | 1817 | 1712
$\mathrm{Growthrate}/\mathrm{(rad/s)}$ | $\times$ | $\times$ | $\times$ | $\times$ | 409.2 | $\times$ | $\times$

The list of model parameters used in the above low-order calculation are given
in Table 1. The length of the plenum is varied while keeping all other
parameter unchanged to investigate its effect on the instability mode
characteristics. Note that the length used for the LES mesh is 0.2 m (see Fig.
1). The modelling result for the frequency and growth rate of the
corresponding mode versus the plenum length is presented in Table 2. The
frequency of 1900 Hz for the plenum length equals 0.2 m agrees very well with
LES result of 1910 Hz. Since the FTF is only reliable at the excited mode
frequency, only the growth rate at 1900 Hz is presented. The predicted growth
rate is positive implying that the system is unstable, which is again
consistent with the LES. It can be also seen that the predicted frequency
monotonically decreases with the increasing plenum length, showing a strong
longitudinal dependency. This provides an explanation for the over-predicted
instability frequency in the LES because an effectively shorter length is
considered without including the bottom conical inlet volume. The height of
this cone is about 60 mm, which would give an effective length of 20 mm if
converted to a cylindrical volume with the same circular area. Consequently,
the corresponding plenum length $L_{0}=0.22$ m gives a frequency of 1817 Hz,
which is close to the measured value around 1800 Hz. This observation suggests
that one should be extremely careful while simplifying the burner geometry for
either LES or low-order modelling as small details may have significant impact
on the system instability behaviour.

Fig. 10: Pressure mode shape at the mean radius of the model combustor when
$L_{0}=0.2$ m. Fig. 11: Axial pressure mode shape at the mean radius of the
model combustor when $L_{0}=0.2$ m and $\theta=0^{\circ}$.

The predicted mode shape by the low-order model is shown in Fig. 10. The
overall distribution also matches well with that obtained by the LES (see Fig.
8). Figure 11 depicts the axial pressure mode shape at the mean radius of the
model combustor when $L_{0}=0.2$ m and $\theta=0^{\circ}$. It is seen that one
wavelength exists in both the plenum and premixer, while the pressure
amplitude in the premixer is larger than that in the plenum. Although the
combustion chamber occupies the largest length scale, the axial pressure
variation there is the smallest due to the higher temperature and hence larger
speed of sound and wavelength. A notable difference between the mode shapes
predicted by LES and low-order model is that the LES gives a higher amplitude
in the plenum compared to that in the combustion chamber, whereas it is the
opposite in the low-order modelling. It is speculated that this is caused by
the simplified outlet boundary geometry condition, i.e. a flat open outlet
shown in Fig. 9, whereas in the present LES and previous Helmholtz simulation
[12] the outer wall is much higher than the inner wall (see Fig. 8). The
former is fully reflective while the latter is not. As a result, the amplitude
in the combustion chamber given by the low-order model is relatively higher
than that obtained from LES and Helmholtz solver. Unfortunately, experimental
measurement of the pressure amplitudes in the plenum and combustion chamber is
unavailable for validation. Nevertheless, the LES and low-order modelling
results presented here provide useful information for the pressure mode
structures inside the combustor, which complements the experimental
observations with useful physical insights. Future works will be directed
towards understanding and predicting the various modal dynamics and explore
the effects of the downstream boundary conditions.

##  6 Conclusions

A comprehensive numerical study including LES and low-order modelling of the
azimuthal mode thermoacoustic instability in an annular combustor with non-
swirling flows is presented. Two LES flamelet models, with and without
considering the heat losses, are used to investigate non-adiabatic effects on
the azimuthal instability. The azimuthal modes are captured by both models and
only marginal differences are observed for the frequency and amplitude. The
ability of flamelet-based LES approach to capture self-excited azimuthal modes
is demonstrated for the first time in this work. Although the amplitudes are
predicted well, the computed frequencies are higher than the measured value
(about 1800 Hz) by 8$-$12%, which is due to the simplified upstream plenum
geometry. The subsequent low-order modelling results show that when the
effectively plenum length is matched the frequency agrees very well with the
measurement.

Analysis of the pressure data from experiments and LES indicates an anti-
clockwise spinning mode with the spin ratio close to positive unity.
Comparison of phase-averaged heat release rate fluctuations shows a good
agreement for this spinning mode. DMD is applied to further investigate the
modal behaviours at the azimuthal mode frequency. A non-zero combustor outlet
mode structure is observed for the pressure in the LES, which suggests that
the downstream geometry may play an important role in establishing the
azimuthal modes in this combustor. The LES observation is further supported by
the low-order modelling results showing very similar mode shapes in the plenum
and premixers, but different in the combustion chamber due to a simplified
outlet boundary. The higher amplitude in the combustion chamber indicates that
a fully reflective open boundary condition may be inadequate and possible
improvements will be explored in future works.

## Appendix A. Assessment of mesh sensitivity

To avoid the high computational cost of repeatedly running full annular
combustor simulations, the mesh sensitivity is investigated using a single
sector of the 12 burner configuration [16]. A schematic of the simulation
setup is shown in Fig. A.12. The bulk velocity at the bluff-body exit is
$U_{b}=27$ m/s and the equivalence ratio is 0.9 for the premixed ethylene/air
mixture at 300 K and 1 atm. Two meshes having 1.2 and 2.1 million tetrahedral
cells are tested here and the typical grid spacing near the shear layer is
about 0.6 and 0.3 mm respectively.

Fig. A.12: Illustration of the single sector simulation with upstream velocity
forcing. Left schematic adapted from [38], copyright 2016, with permissions
from ASME.

First, the unforced reacting flow is simulated for the 1.2M and 2.1M meshes.
After the flame is fully established, statistics are collected for about 4
flow-through-time, which is calculated based on the combustion chamber length
of 80 mm and the bulk velocity. Figure A.13 compares the computed mean axial
velocity and temperature radial profiles at axial distances of $x/D=0.5$, 1
and 2. The bluff-body diameter is $D=13$ mm. As can be seen, the two
simulations give very similar results for the mean axial velocity field at all
three streamwise positions except for some considerable differences in the
bottom of the recirculation zone. This is probably due to the relatively short
statistics sampling time for the 2.1M case as indicated by the asymmetry. The
difference in the temperature fields is almost negligible suggesting a good
grid convergence.

Fig. A.13: Mean axial velocity and temperature comparisons for the unforced
case.

Further grid sensitivity is assessed by performing LES with forced inlet
velocity (See Fig. A.12) to compare the flame transfer functions given by the
two meshes. As the azimuthal modes are of interests here, a fixed forcing
frequency of 1800 Hz corresponding to the first azimuthal mode based on the
perimeter of annular chamber is considered. Also, as the amplitude of pressure
oscillation is below 1% of the ambient pressure [12], the two relatively low
forcing amplitudes, 5% and 10% of the mean flow velocity respectively, are
used for the forcing.

Fig. A.14: Time series of the fluctuating velocity at the bluff-body exit
(black) and the heat release response (red) to the excitation for the 1M
(solid) and 2M (dashed) meshes.

Figure A.14 presents the time series (after transition) of the normalised
fluctuating velocity at the bluff-body exit and heat release response to the
excitation for the two meshes (see Fig. A.12 for the formulae used). The
overall fluctuating behaviour of the two signals are quite similar despite
some slight differences in the heat release amplitude and the time delay
(within 10%). To show this in clearer manner, the gain and phase of the
transfer functions are computed as shown in Fig. A.15. Again, similar values
are obtained from the two meshes for both the gain and phase. From 5% to 10%
forcing amplitude, the gain slightly decreases and the phase increases, and
this trend is also captured in both the 1.2M and 2.1M cases.

The above grid sensitivity study on the single sector suggests that the
smaller 1.2M mesh is sufficient to capture the reactive flow fields and also
the flame response to external forcing excitations. Thus, this grid
distribution is applied to all 12 sectors for the annular combustion LES shown
in the main text of this paper.

Fig. A.15: Time series of the fluctuating velocity at the bluff-body exit
(black) and the heat release response (red) to the excitation for the 1M
(solid) and 2M (dashed) meshes.

## Acknowledgements

ZXC and NS acknowledge the support of Mitsubishi Heavy Industries, Takasago,
Japan. This work used the ARCHER UK National Supercomputing Service
(http://www.archer.ac.uk) with the computational time provided by the UK
Consortium for Turbulent Reacting Flow.

## References

  * [1] T.C. Lieuwen, V. Yang (Eds.), Combustion instabilities in gas turbine engines: Operational experience, fundamental mechanisms, and modeling, AIAA Inc., Reston, VA, 2006.
  * [2] T. Poinsot, Prediction and control of combustion instabilities in real engines, Proc. Combust. Inst. 36 (2017) 1–28.
  * [3] J. O’Connor, V. Acharya, T. Lieuwen, Transverse combustion instabilities: Acoustic , fluid mechanic , and flame processes, Prog. Energy Combust. Sci. 49 (2015) 1–39.
  * [4] J.R. Dawson, N.A. Worth, Flame dynamics and unsteady heat release rate of self-excited azimuthal modes in an annular combustor, Combust. Flame 161 (2014) 2565–2578.
  * [5] N.A. Worth, J.R. Dawson, Modal dynamics of self-excited azimuthal instabilities in an annular combustion chamber, Combust. Flame 160 (2013) 2476–2489.
  * [6] J.F. Bourgouin, D. Durox, J.P. Moeck, T. Schuller, S. Candel, A new pattern of instability observed in an annular combustor: The slanted mode, Proc. Combust. Inst. 35 (2015) 3237–3244.
  * [7] L. Crocco, Aspects of combustion instability in liquid propellant rocket motors. Part I: Fundamentals. Low frequency instability with monopropellants, J. Am. Rocket Soc. 21 (1951) 163–178.
  * [8] L. Crocco, Aspects of combustion instability in liquid propellant rocket motors. Part II: Low frequency instability with bipropellants. High frequency instability, J. Am. Rocket Soc. 22 (1952) 7–16.
  * [9] M. Bauerheim, F. Nicoud, T. Poinsot, Progress in analytical methods to predict and control azimuthal combustion instability modes in annular chambers, Phys. Fluids 021303.  doi:10.1063/1.4940039.
  * [10] A.P. Dowling, S.R. Stow, Acoustic analysis of gas turbine combustors, J. Propul. Power 19 (2003) 751–764.
  * [11] X. Han, J. Li, A.S. Morgans, Prediction of combustion instability limit cycle oscillations by combining flame describing function simulations with a thermoacoustic network model, Combust. Flame 162 (2015) 3632–3647.
  * [12] M. Bauerheim, G. Staffelbach, N.A. Worth, J.R. Dawson, L.Y.M. Gicquel, T. Poinsot, Sensitivity of LES-based harmonic flame response model for turbulent swirled flames and impact on the stability of azimuthal modes, Proc. Combust. Inst. 35 (2015) 3355–3363.
  * [13] P. Wolf, G. Staffelbach, L.Y.M. Gicquel, J.-D. Müller, T. Poinsot, Acoustic and Large Eddy Simulation studies of azimuthal modes in annular combustion chambers, Combust. Flame 159 (2012) 3398–3413.
  * [14] N. Zettervall, N.A. Worth, M. Mazur, J.R. Dawson, C. Fureby, Large eddy simulation of CH4-air and C2H4-air combustion in a model annular gas turbine combustor, Proc. Combust. Inst. 37 (2019) 5223–5231.  doi:10.1016/j.proci.2018.06.021.
  * [15] J.-F. Bourgouin, D. Durox, J.P. Moeck, T. Schuller, S. Candel, Self-sustained instabilities in an annular combustor coupled by azimuthal and longitudinal acoustic modes, in: Proc. ASME Turbo Expo, no. GT2013-95010, Texas, USA, 2013\.  arXiv:https://asmedigitalcollection.asme.org/GT/proceedings-pdf/GT2013/55119/V01BT04A007/4223970/v01bt04a007-gt2013-95010.pdf, doi:10.1115/GT2013-95010.
  * [16] M. Mazur, H.T Nygard, J.R. Dawson, N.A. Worth, Characteristics of self-excited spinning azimuthal modes in an annular combustor with turbulent premixed bluff-body flames, Proc. Combust. Inst. 37 (2019) 5129–5136.
  * [17] C. Kraus, L. Selle, T. Poinsot, Coupling heat transfer and large eddy simulation for combustion instability prediction in a swirl burner, Combust. Flame 191 (2018) 239–251.
  * [18] P.W. Agostinelli, D. Laera, I. Boxx, L. Gicquel, T. Poinsot, Impact of wall heat transfer in Large Eddy Simulation of flame dynamics in a swirled combustion chamber, Combust. Flame 234 (2021) 111728.  doi:10.1016/j.combustflame.2021.111728.
URL https://doi.org/10.1016/j.combustflame.2021.111728

  * [19] J.C. Massey, Z.X. Chen, N. Swaminathan, Modelling heat loss effects in the Large Eddy Simulation of a lean swirl-stabilised flame, Flow, Turbul. Combust. 106 (2021) 1355–1378.
  * [20] Z.X. Chen, N. Swaminathan, M. Stöhr, W. Meier, Interaction between self-excited oscillations and fuel-air mixing in a dual swirl combustor, Proc. Combust. Inst. 37 (2019) 2325–2333.
  * [21] Z.X. Chen, I. Langella, N. Swaminathan, M. Stöhr, W. Meier, H. Kolla, Large Eddy Simulation of a dual swirl gas turbine combustor: Flame/flow structures and stabilisation under thermoacoustically stable and unstable conditions, Combust. Flame 203 (2019) 279–300.
  * [22] J. Su, A. Garmory, J. Carrotte, The aerodynamic response of fuel injector passages to incident acoustic waves, in: Proc. ASME Turbo Expo, no. GT2015-43248, Montréal, Canada, 2015.
  * [23] J. Smagorinsky, General circulation experiments with the primitive equations, Mon. Weather Rev. 91 (3) (1963) 99–164.
  * [24] Z. Chen, S. Ruan, N. Swaminathan, Large Eddy Simulation of flame edge evolution in a spark-ignited methane-air jet, Proc. Combust. Inst. 36 (2017) 1645–1652.
  * [25] H. Pitsch, H. Steiner, Large-eddy simulation of a turbulent piloted methane/air diffusion flame (sandia flame D), Phys. Fluids 12 (2000) 2541.
  * [26] I. Langella, Z.X. Chen, N. Swaminathan, S.K. Sadasivuni, Large-Eddy Simulation of reacting flows in industrial gas turbine combustor, J. Propul. Power 34 (2018) 1269–1284.
  * [27] H. Wang, A. Laskin.  A comprehensive kinetic model of ethylene and acetylene oxidation at high temperatures [online] (1998).
  * [28] C. Xu, A.A. Konnov, Validation and analysis of detailed kinetic models for ethylene combustion, Energy 43 (2012) 19–29.  doi:https://doi.org/10.1016/j.energy.2011.11.006.
  * [29] R. Balachandran, B.O. Ayoola, C.F. Kaminski, A.P. Dowling, E. Mastorakos, Experimental investigation of the nonlinear response of turbulent premixed flames to imposed inlet velocity oscillations, Combust. Flame 143 (2005) 37–55.
  * [30] X. Han, A.S. Morgans, Simulation of the flame describing function of a turbulent premixed flame using an open-source LES solver, Combust. Flame 162 (2015) 1778–1792.
  * [31] S. Ruan, T.D. Dunstan, N. Swaminathan, R. Balachandran, Computation of forced premixed flames dynamics, Combust. Sci. Technol. 188 (2016) 1115–1135.
  * [32] G. Ghirardo, M.P. Juniper, Azimuthal instabilities in annular combustors: standing and spinning modes, Proc. R. Soc. A 469 (2157) (2013) 20130232.
  * [33] G. Ghirardo, F. Gant, Background noise pushes azimuthal instabilities away from spinning states, arXiv preprint arXiv:1904.00213.
  * [34] P.J. Schmid, Dynamic mode decomposition of numerical and experimental data, J. Fluid Mech. 656 (2010) 5–28.  doi:10.1017/S0022112010001217.
  * [35] C. Huang, W.E. Anderson, M.E. Harvazinski, V. Sankaran, Analysis of self-excited combustion instabilities using decomposition techniques, AIAA J. 54 (9).
  * [36] G. Zhang, X. Wang, L. Li, X. Sun, Effects of perforated liners on controlling combustion instabilities in annular combustors, AIAA J. 58 (7) (2020) 3100–3114.  doi:10.2514/1.J059099.
  * [37] S Akamatsu, A P Dowling, Three dimensional thermoacoustic oscillation in a premix combustor, in: ASME Turbo Expo 2001: Power for Land, Sea, and Air, American Society of Mechanical Engineers, 2001, pp. V002T02A001–V002T02A001.  doi:10.1115/2001-GT-0034.
  * [38] G. A. Mensah, G. Campa, J. P. Moeck, Efficient Computation of Thermoacoustic Modes in Industrial Annular Combustion Chambers Based on Bloch-Wave Theory, J. Eng. Gas Turbines Power 138 (8) (2016) 1–7.  doi:10.1115/1.4032335.