Azimuthal flame response and symmetry breaking in a forced annular combustor

In the current study azimuthal forcing of an annular combustor with swirling flames has been performed to present for the first time the Heat Release Rate (HRR) response to all possible pressure fields of the first azimuthal mode up to a finite amplitude limit. The response is first quantified through the conventional Flame Describing Function (FDF) framework, showing a difference in response which depends on whether the acoustic field rotates anti-clockwise or clockwise, albeit with some scatter. Additionally and somewhat surprisingly, a finite HRR response is observed for flames exactly in the pressure node. An Azimuthal FDF is introduced, based on the decomposition of the spatial HRR response through Bloch theory, to better highlight the difference in HRR response to the anti-clockwise and clockwise components of the acoustic field and reduce scatter. A clear difference in response is observed, with a significantly higher response to the anti-clockwise forcing component compared to the clockwise component, independent of the prescribed pressure mode. The difference is attributed to the systematic symmetry breaking introduced by using an annular enclosure of finite curvature and width with swirling flames. It is argued that the finite curvature and width of the geometry and the swirl need to be both present to observe this effect. The difference in response results in a difference between the nature angle of the HRR mode and that of the acoustic field, explaining the relatively large HRR response observed for flames in the pressure node. The Azimuthal FDF describe all of these phenomena well, and is therefore considered more suitable than the conventional FDF to characterise the response in an annular combustor. © 2021 The Author(s). Published by Elsevier Inc. on behalf of The Combustion Institute. This is an open access article under the CC BY license ( http://creativecommons.org/licenses/by/4.0/ )


Introduction
Thermoacoustic instabilities are a well known development issue for gas turbine engines [1] . These occur due to the interaction between heat release rate and pressure oscillations, and are prevalent when operating over very wide ranges of conditions or when varying fuels; both of which may be advantageous for the control of emissions in modern engines. Thus, in order to realise the benefits of increased fuel and operational flexibility, accurate prediction of thermoacoustic stability is important.
One promising approach to predict the thermoacoustic stability of a system is to use either low-order acoustic network models [2][3][4][5][6] or Helmholtz solvers [7][8][9] . These employ Flame Transfer or Flame Describing Functions (FTF/FDF) in order to couple the nonlinear flame response to the acoustics of the system, and are usually obtained experimentally [10,11] , or numerically [12,13] . While can also vary slowly with time, resulting in distinct modal dynamics [14,15,20,21] . Furthermore, the presence of multiple flames around an annular chamber means that neighbouring flames are free to interact with each other [22] , adding further complexity to the response.
At present, FDFs are commonly measured, simulated or defined based on isolated single flame setups subjected to acoustic oscillations in the axial, or bulk flow, direction [10][11][12][13]23] . Such functions are therefore suitable for predicting the stability of single sector rigs in order to demonstrate the methodology [24][25][26] , and may also be suitable for describing the stability of more practical can-annular configurations [27] . However, transfer functions obtained from single sector measurements omit by definition features associated with the more complex response found in annular configurations, and include no functional dependence on χ or θ 0 , which is a potential source of uncertainty in such stability predictions. While a handful of studies have used either simplified transfer function models, or functions obtained from single sector measurements to predict the stability features of an annular system [9,[28][29][30][31][32][33][34] , the use of single flame transfer functions means any possible effects due to adjacent flames interacting with each other are eliminated. Additionally, any effects on the transfer function from the bulk swirl, due to the geometry, combined with the azimuthal mode are not included. It therefore remains an open question whether transfer functions obtained in single flame setups are more generally applicable for stability prediction in annular configurations, as the functional dependence of these on χ or θ 0 is as yet unknown.
In order to better understand the flame response to azimuthal oscillations, a number of studies have examined single flames confined within enclosures which are long in the transverse direction (orthogonal to the bulk flow) [35][36][37][38][39][40] . Such arrangements are useful as they allow the effect of transverse acoustic excitation to be studied across a wide range of frequencies, while removing the complexity associated with flame-flame interactions. More recently, transverse oscillations have also been investigated in linear arrays of multiple flames [41] , permitting the effect of flame-flame interactions to be studied in addition to the effect of transverse oscillations on the flame dynamics. However, all systems with linear enclosure geometry omit effects associated with finite wall curvature combined with a bulk swirl in the annulus induced by individual swirling flames, prohibiting the investigation of these symmetry breaking effects.
A number of recent theoretical studies have explored symmetry breaking effects, through the derivation of dynamical equations for the azimuthal mode state space variables [42,43] . The degenerate eigenvalue has been shown to split when the reflectional symmetry is broken [44] , leading to two different growth rates [43,45] , but these effects are not related to changes in the heat release rate. While the describing function in these dynamical formulations is allowed to have a dependence on the nature angle χ, this has not been widely explored, as there is as yet no systematic evidence linking the nature angle and heat release rate. Therefore, while it is widely known that azimuthal modes can contain a rich range of model responses, and there exists a way to include these in dynamical equations governing the system stability, there is as yet little guidance on how to link the dependence of the response function to the azimuthal mode nature.
In order to address this the current study is aimed at capturing for the first time the full dependence of flame describing functions in an annular geometry on the nature, orientation, and amplitude of azimuthal modes. In order to exert control over these parameters, the excitation strategy of [46] is applied in an annular combustor. This approach has been previously demonstrated for the excitation of standing [19] and single amplitude spinning modes [47] , but the full range of spinning and mixed standing and spinning modes have not been explored yet. The present study realises a 40 fold increase of the number of parameter combinations compared to previous work, which allows a detailed examination of the functional dependence of the FDF on nature angle, orientation angle and amplitude leading to new and somewhat unexpected results. Different spatial locations within a transverse standing mode [4 8,4 9] and a transverse travelling acoustic mode [38] have been studied in single flame setups, and it has been shown that small asymmetries in single flames can cause a difference in response when subjected to transverse acoustic velocities [50][51][52] . However, this study is the first to examine the effect of the systematic symmetry breaking associated with the presence of swirling flames in an annular combustion chamber. The analysis is performed in a way that mini-mizes the effects of small flame to flame differences, shifting the focus from the effect on a single flame to the effect on the global response.
It is also worth noting that these measurements also permit examination of flame to flame differences, by comparing nominally identical forced states on 12 nominally identical flames positioned around the annular chamber. Finite physical differences between individual flames are another source of symmetry breaking, which has been shown to affect system stability [53] . To analyse this, a novel implementation of the Bloch formalism [8,54] is introduced, enabling the calculation of an average flame (which is independent of flame to flame differences), and the decomposition of the first azimuthal mode into the individual ACW and CW components. In this manner the FDF of a single, average flame can be defined. Additionally the decomposition is used to define a so called Azimuthal Flame Describing Function ( F DF ± ) for each component spinning ACW ( -) and CW ( + ) of the mode, which is shown to be the most descriptive representation of the heat release rate response in the current annular configuration. The Azimuthal FDFs are also used to provide a theoretical, nature angle dependent HRR description, suitable for use in the recently derived dynamical equations [42,43] .
The rest of the paper is organised as follows. First the experimental setup and post processing is described in detail in Section 2 , starting with the setup ( Section 2.1 ), pressure calculation methods ( Section 2.2 ), conventional Flame Describing Function definition ( Section 2.3 ), how to prescribe the pressure modes ( Section 2.4 ), phase averaging ( Section 2.5 ), before finally describing rotational averaging and the corresponding Bloch waves formalism ( Section 2.6 ) and introducing the concept of the HRR mode ( Section 2.7 ). The measured FDFs of the individual flames are presented in various forms in Sections 3.1 -3.3 , followed by the two Bloch kernels in Section 3.4 . The new Azimuthal FDFs are presented in Section 3.5 , before the implications on the HRR modes are discussed and quantified in Section 3.6 . Lastly the conclusions are given in Section 4 .

Geometry and data acquisition
The combustor used in this study is the annular combustion chamber used in [47] in the N = 12 flame configuration, as shown in Fig. 1 . The combustor is described in detail in earlier studies, but the main features and dimensions will be repeated here. The combustion chamber consists of two concentric cylindrical walls of diameter d outer = 212 mm and d inner = 129 mm . To be comparable to the self-excited configuration, and help with optical access from above, the inner wall is shorter than the outer wall, with lengths of 127 mm and 287 mm respectively. Both walls are mounted at the dump plane, which is defined as the plane separating the injector tubes and the combustion chamber as shown in Fig. 1 c. A speaker array, consisting of 8 equidistantly spaced, 68 mm long standoffs, is mounted with the speakers approximately 104 mm above the dump plane, ensuring the speakers do not interact directly with the flames through velocity oscillations induced at the exits of the standoff tubes. Four of the standoffs are equipped with Adastra HD60 horn drivers, which are used to prescribe the forced state on the system. The remaining standoffs are left blocked at the outer end of the array, marginally improving the control of the resulting forced state. The horn drivers are driven by an Aim-TTI TGA 1244 signal generator, amplified by a pair of QTX PRO 10 0 0 power amplifiers.
Each of the 12 equidistantly spaced injectors have a diameter of 19 mm , with a central bluff body of diameter 13 mm and half angle 45 • . The bluff body is mounted to a rod of diameter 5 mm used for centering and acts as a mounting point for a swirler. The swirler is mounted 10 mm below the dump plane, measured from the trailing edge of the swirler vanes, and induces an ACW swirl when viewed from downstream. The resulting swirl number 10 mm downstream of the dump plane has been measured to be 0 . 65 for an unconfined configuration. The injectors are 150 mm long and they are fed by a plenum with flow straighteners and a hemispherical body for equal flow distribution. The plenum is fed by two impinging jets of premixed air and ethylene, at an equivalence ratio of = 0 . 7 . The flow rate and equivalence ratio is controlled by three Alicat mass flow controllers, and the flow rate is set to give a bulk velocity of U bulk = 18 m / s at the dump plane. This results in a thermoacoustically stable operating condition for the studied configuration, easing the control of the forced azimuthal mode.
Six of the injector tubes are instrumented with two differential pressure transducers of the type Kulite XCS-093-0.35D flush mounted with the injector wall, spaced by 65 mm in the flow direction. The signals are amplified by a pair of Fylde FE-579-TA bridge amplifiers and digitized by a set of NI-9234 DAQ modules, operating at 51 . 2 kHz . Heat release rate data is obtained by a Phantom V2012 high speed camera equipped with a LaVision Intensified Relay Optics (IRO) unit. The IRO is equipped with a Cerco 2178 UV lens with a narrow band pass filter, centered at 310 nm with a full width half maximum of 10 nm . This captures the light intensity with a wavelength corresponding to de-exciting OH * radicals, which has been shown to be proportional to the heat release rate for perfectly premixed combustion [55] . The imaging system operates at 10 kHz , which is sufficient time resolution for the forced states of interest, and the IRO gate time is 80 μs. The sample rate of the system is not a multiple of the frequency of the forced state by design, improving the number of phase instances captured by the camera. The lengthened trigger signal of the IRO unit is acquired on the same system as the pressure transducers, enabling synchronization of the pressure and heat release rate signals. A spatial resolution of 2.5 pixels per millimeter is achieved and a minimum of 20,0 0 0 images are taken for each forced state. A total of 5,250,0 0 0 forced flame images were taken as part of the current study.

Mode reconstruction
Due to the long aspect ratio of the injector tubes, and small diameter relative to the wavelength of the acoustic mode, the acoustic pressure fluctuations in the tube are two counter propagating oscillations at all azimuthal locations. Amplitudes range from 1224 Pa for the High forcing level, through 661 Pa ( Middle ) and 340 Pa ( Low ) to the noise level ( < 10 Pa ) for the thermoacoustically stable Unforced case. Due to the sharp peak at the prescribed forcing frequency, the frequency of oscillation will be approximated to be constant. plane waves, where denotes the real part of its argument { . . . } . The amplitude of the two counter propagating components are given by B + and B − , corresponding to the component propagating in the downstream and upstream direction respectively. x is the vertical position in the tube as shown in Fig. 1 . The wavenumbers k ± are given by [56] where c is the speed of sound in the medium and Ma is the Mach number. The left hand side of Eq. (1) is measured at two discrete locations, x upper = −44 mm and x lower = −109 mm , by the pressure transducers in the instrumented injectors. The analytical signal of p d (t) is then obtained through the Hilbert transform H of the measured signals, The acoustic velocity perturbations corresponding to the pressure in Eq. (1) can be derived from the fluctuating momentum equation, and are given by [56] The expression for the velocity in Eq. (3) is calculated from the pressure measurements at the two microphone locations by solving Eq. (1) for B ± e i ωt . In general, the solution of Eqs. (1) and (3) is performed for all frequencies, but as a computational simplification the prescribed forcing frequency is assumed to be the only frequency of significance. This assumption is based on the amplitude spectrum shown for a few amplitudes of the forced states and an unforced case in Fig. 2 . The amplitude at the forcing frequency of f = 1650 Hz is several orders of magnitude above the background noise level, and the first harmonic is about two orders of magnitude lower for the presented cases. Combined with the use of the same equivalence ratio and inlet velocity for all forced states, the assumption of a single frequency is reasonable. For a constant cross section, propagating the pressure and acoustic velocity fluctuations to alternate locations is performed by inserting the location x into Eqs. (1) and (3) . In the case of smooth area changes, A 1 → A 2 , the pressure and mass are conserved across the jump in the zero Mach number limit [2] p d The swirler will introduce an initial area decrease, followed by an increase back to the original area. This temporary area change, and other potential effects of the swirler on the axial acoustic mode, are included through the use of the experimentally measured scattering matrix [57,58] . The area change introduced by the bluff body is estimated as 10 separate step changes in the diameter, each with the same length and volume as the corresponding section of the injector and bluff body system. After propagating the pressure oscillations to the dump plane of the combustor, the injector tubes terminate in an annular geometry with an azimuthal pressure mode. The acoustic mode in the injector tubes couples primarily to the acoustic pressure in the chamber and not to the acoustic azimuthal velocity [59,60] . Since all the injectors are nominally identical, and the potential dependence of the impedance on the mode is low [59] , it is safe to assume the impedance linking the pressure in the injectors and dump plane is the same for all the injectors. Neglecting the specific value of the pressure coupling impedance, the pressure in the combustion chamber can be described by the following pressure ansatz describing the first azimuthal mode [18] p ( θ , t ) = A cos ( θ − θ 0 ) cos (χ ) cos ( ωt + ϕ ) The derivation based on quaternions and detailed procedure for fitting the measured pressure at the dump plane can be found in [18] . In the current study a least squares solution based on the six propagated pressures are used to reduce the effect of random fluctuations in the pressure signal. The amplitude of the azimuthal pressure mode is denoted A in Eq. (5) , and χ is the nature angle of the mode. The nature angle describes whether the mode is purely spinning ( χ = π / 4 for ACW, χ = −π / 4 for CW), purely standing ( χ = 0 ), or a mix of the two for intermediate values. It is equivalent to the Spin Ratio (SR) [15] used in previous studies [14,47] on the same geometry through the relation SR = tan (χ ) . The orientation angle θ 0 gives the location of the pressure anti-node, and ϕ is a temporal phase. The pressure fluctuations in Eq.
The characteristic features of the pressure mode for a spinning and a standing mode are shown for two example forced states in Fig. 3 . The spinning mode in the top row show all the pressure locations have similar amplitude, and the phase difference corresponds to the time it takes the wave to cover the spatial distance between the pressure transducers. The standing mode is in con-trast observed to have two similar amplitude pressure signals of opposite phase, with a third signal of negligible amplitude. This is caused by the characteristic alternating pattern of nodes and antinodes, equidistantly distributed every azimuthal angle π / 2 for the first azimuthal mode. For the range of nature angles χ between the standing and the spinning modes, the mode can be considered a mixed mode, with a standing component and a spinning component. The concept of the orientation angle θ 0 only makes sense in the case of mixed or standing modes where A st = 0 .

Flame describing functions (FDFs)
The interaction between a single flame and acoustics are often described by the conventional Flame Describing Function ( FDF ) [61] of a single flame, where the spatially integrated HRR of the flame is considered. In the current configuration there are N = 12 different flames inside the combustion chamber. To extract the response most equivalent to a single flame setup, the response within one flame sector , as defined in Fig. 1 b, is considered. The spatially integrated HRR of each flame sector, denoted q j for sector number j, is obtained by summing over the pixels within the sector. This removes the spatial dependence within the jth flame sector, quantifying the HRR of a single flame as a scalar. The conventional FDF for the jth individual flame in the annular combustor from the injector located at θ = θ j is then defined as [61] F DF j ω, ˆ where ˆ u axial , j is the axial velocity perturbation from Eq. (3) in the jth injector at a reference axial location, here chosen to be the dump plane. ˆ q j is the complex Fourier amplitude of the heat release rate oscillations at the prescribed angular excitation frequency ω = 2 π f , and the mean heat release rate of the jth flame is denoted q j . In recent work [43,53,62] similar approaches have been used to describe the interaction between heat release rate and pressure in annular configurations in the governing equations. Usually the describing function is obtained over a range of frequencies f for longitudinal forcing setups. However, in this study only the response of the first azimuthal mode is of interest, meaning the frequency will be fixed by the geometry and temperature of the combustion chamber. The response is still defined as a describing function [63,Section 2.3] , and the naming is in line with the first use of describing functions in thermoacoustic instability studies by Dowling [3,Eq. (3.8)] . The study of a single frequency is in contrast to one previous study on a similar configuration [64] , where a large range of frequencies were used, and the nature of the excited modes was not explicitly controlled, making the results more difficult to interpret. Instead in the present study, the describing function at a fixed frequency will be explored for different combinations of the state space variables A , χ and θ 0 in Eq. (5) for six different flames directly from experimental measurements, and from 12 flames by interpolation of the pressure wave oscillations. This is the first time the describing function is investigated for 12 nominally identical swirling flames arranged in an annulus, enabling assessment of the variability between the injectors. Therefore, the flame to flame differences can be assessed in detail.

Prescribing forced states
Previous experimental studies featuring self-excited azimuthal instabilities in annular combustors have not been able to systematically explore the full parameter space provided by Eq. (5) due to constantly changing state space parameters [14] and statistical preference for certain parameter combinations [19] . Therefore, Fig. 4. All the 123 forced states on the Poincaré sphere, with 3D representation ( middle left ), projection on the plane χ = 0 ( middle right ) and projection on the θ = 0 plane ( right ). The relation between the parameters of the forced state and the coordinate system is illustrated on the left . The orientation angle θ 0 is restricted to the interval θ ∈ [ 0 , π ) . a forcing array is used to focus on the HRR for a fixed system state, expanding the limited combination of states that exist in self-excited systems. Forcing arrays have been used previously for limited combinations of state space variables in annular swirling combustors [19,47,64] , and in the current study the application of the technique is improved to obtain well controlled forced states for the full parameter space up to the physical amplitude limit of the forcing array. Similarly to the previous study in [47] , two pairs of horn drivers are used to impose the forced state. The two horn drivers in each pair are located on diametrically opposite sides of the combustion chamber, and are driven 180 • out of phase. The two pairs are separated by an azimuthal angle of 90 • .
The strategy for imposing forced states follows the simple idea of matching the amplitude of both speaker pairs, each setting up a separate standing mode, and then adjusting the phasing between the two pairs to obtain the desired nature angle χ and orientation angle θ 0 . The forcing array is applied to a thermoacoustically stable operating point to ensure good control of the mode. While the mode is thermoacoustically stable, the preference for certain nature angles observed in self-excited studies is still prevalent due to the reacting flow [19] . Additionally, the combustor can only be run for a limited time due to the temperature limits of the forcing array, meaning total thermal equilibrium cannot be achieved.
To use a fixed frequency of f = 1650 Hz the combustor is always ignited at similar reference outer wall temperature, and the data acquisition is started approximately 20 s later, ensuring the temperatures in the chamber are highly repeatable, and well-tuned to the forcing frequency. The frequency is determined by iteration over different frequencies at constant horn driver power, choosing a frequency close to the maximum response at the upper pressure transducer location while exhibiting good control of the mode. The increasing temperature of the combustion chamber after ignition introduces a slight drift in the mode from ignition to the time of data acquisition. However, the frequency is chosen such that the drift is negligible during the measurement period. A live pressure mode reconstruction is used to constantly monitor the mode, making it possible to make adjustments to the forcing array until the measurements start. In general, standing forced states χ ≈ 0 are very easy to obtain experimentally, while the spinning forced states χ ≈ ±π / 4 are relatively hard to obtain experimentally. This difference is most likely caused by combination of the preference for the standing modes due to noise [42] , and a requirement for a delicate amplitude balance between the speaker pairs for the spinning modes.
The final 123 prescribed forced states are indicated on the Poincaré sphere [18] in Fig. 4 , and are summarised briefly in Table 1 . The first azimuthal mode is the only one considered, and excited, in the current work. The maximum acoustic pressure level is of the order of 1% of the operating atmospheric pressure, representative of industrial applications [65] . The number of Table 1 Overview of the range for the pressure mode parameters for the studied forced states. Close to the full parameter space is studied up to a finite amplitude limit imposed by the speaker array.

Parameter
Range forced states is more than 40 times higher than previous studies in similar configurations, and higher magnitude nature angles are achieved [19,47] in addition to a range of different amplitudes. Standing modes, as presented in the lower part of Fig. 3 , live in the horizontal plane cutting the center of the sphere in Fig. 4 . The angle in the horizontal plane corresponds to the orientation angle θ 0 of the pressure anti-node location. Conversely the spinning modes, for example ACW spinning mode shown in the uppert plot of Fig. 3 , live along the vertical axis in Fig. 4 . Everything in between are called mixed modes, and contain a mix of a standing component and a spinning component. Since the orientation angles θ 0 and θ 0 + π of the standing mode component are equivalent, the range 0 − π is the main focus. The standing and mixed forced states are repeated for each π / 4 increment in orientation angle in this interval, ensuring each injector is subjected to different points in the standing component of Eq. (5) . The final number of forced states comes from a desire to have at least 5 different amplitudes for the 4 unique nodal line positions for all the standing and mixed mode nature angles, as well as different amplitudes for the spinning modes. This is the first time this level of control and systematic exploration of the state space has been achieved for annular combustors with swirl.

Phase averaging
Turbulent fluctuations play a significant role in the instantaneous heat release rate q , complicating the interpretation of the heat release rate fluctuation structures. In the following, instead of considering a real valued, time dependent function q ( r, θ , t ) , the respective analytic signal q a ( r, θ , t ) will be considered. The analytic signal is defined as q a ( r, H is the Hilbert transform. The analytic heat release rate expression is then given by the spatially dependent mean heat release rate q ( r, θ ) , the phase dependent fluctuating component q a ( r, θ , t ) and the stochastic fluctuating part q s , a ( r, θ , t ) Here the subscript "a" denotes it is the complex valued analytic expression of the fluctuations. The phase dependent fluctuations, which are periodic in time t, can be obtained from phase averaging the response for a sufficiently large number of samples M, with t 0 ∈ [ 0 , T ) where T is the oscillation period. This is due to the stochastic fluctuations having a zero mean for sufficiently many samples The oscillation cycle is divided into 36 equally wide phase bins, resulting in a minimum number of samples of at least M ≈ 550 in each bin, which is considered sufficient. Each sample is then included in the phase average bin closest to the phase provided by the upper pressure transducer ( x = x upper ) at the θ = 0 location, as defined in Fig. 1 .

Rotational averaging and Bloch theory
One defining feature of this forced study is the use of N = 12 injectors arranged equidistantly around the annulus. This makes it a prime candidate to use the concept of Bloch theory [54] , which was recently introduced for simulations of annular combustors by Mensah and Moeck [8] . Mensah and Moeck used this to reduce the computational domain by exploiting the N fold symmetry of a combustor with N sectors to only calculate the response of a single sector. However, in the current study the response of all the flames are measured, and Bloch theory is instead used to find the response of an average flame , which corresponds to the Bloch wave part of the response. This will be done in terms of the so called rotational averaging procedure first introduced in [47] . While the exact same rotational averaging procedure used in [47] is followed here, this will now be interpreted in terms of Bloch theory. It is useful to describe it in such terms, both because the notation of the process can be defined more precisely, and also because of the current paper's focus on the flame to flame differences, which can be defined explicitly using this formulation. It also provides a useful experimental reference for a method which has so far only been applied numerically in the context of thermoacoustic instabilities.

Temporal mean heat release rate
The base idea of the rotation averaging procedure is to average all N = 12 distinct flames together, to create an average flame . In the simplest case this can be used to obtain the temporal mean HRR of the average flame from the temporal mean flames q ( r, θ ) where the N fold rotational symmetry is exploited to add all the distinct flames together to get an average response and no Bloch theory is required. This reduces the spatial dependence to a single flame sector. The temporal mean HRR of the distinct flames q ( r, θ ) and the corresponding temporal mean HRR of the average flame q ( r, θ ) sectors are shown for the chosen operating condition, which is thermoacoustically stable, for an unforced case in Fig. 5 . Overall, the main flame structures are observed to be similar for the two quantities, but there are also some differences in the spatial distribution of the heat release rate from flame to flame. All the flames exhibit both positive and negative differences simultaneously, meaning the net result, shown on the bottom right, are relatively small but should still be considered. Therefore, in the case of examining the first azimuthal wave response of the heat release The net differences ( lower right ), obtained by sector averaging, are presented as colored regions. The structure is observed to be similar for the distinct and average flames, but there are some significant differences in the spatial distribution of the HRR within the individual flames. However, the differences in all the flames are both positive and negative, resulting in a much smaller net difference.
rate fluctuations, the fluctuations should be normalised against the corresponding average flame temporal mean HRR q ( r, θ ) sectors , and not the temporal mean of the distinct flames q ( r, θ ) .

Heat release rate fluctuations
The same fundamental idea of how to obtain an average flame response is also used to introduce the rotation averaging procedure for the fluctuating part of the heat release rate. The base assumption will be that the flames are the same and the only significant azimuthal component in the measured response is the first, degenerate, azimuthal mode. This is reasonable because the prescribed forced state is restricted to the first azimuthal mode. The flames will also be modeled as responding only at the forcing frequency f = ω/ 2 π , because in Fig. 2 the amplitude of the harmonics are at least two orders of magnitude lower for all cases. Utilizing Bloch theory, the analytical phase averaged heat release rate in Eq. (9) is modeled as The functions ψ ±1 ( r, θ ) are the spatially dependent heat release rate mode shapes, which are called Bloch kernels. The corresponding components in Eq. (12) account for the degeneracy of the first azimuthal mode. ψ − is the amplitude of the ACW rotating component, and ψ + is the amplitude of the CW rotating component. The flame to flame differences, that are present in all physical systems, are accounted for in the ˆ q a term. These differences are assumed to be small, made formally explicit by the inclusion of the factor 0 < 1 in the expression. The ˆ q a term in general accounts for the violation of the main assumptions made before Eq. (12) , making it possible to describe an arbitrary response. Eq. (12) is also able to account for flame-flame interactions, assuming that the interaction is the same between the different flames.
The Bloch kernels ψ ±1 are defined on a single flame sector, as defined in Fig. 1 , and are periodic [8] : Inspired by Eq. (11) the rotational average in the ACW (-) and CW (+) directions are introduced as the sum (14) corresponding to rotating the coordinate system an angle ∓2 π /N for each step in time 2 π / ( ωN ) . This definition is exactly the same as originally proposed in [47] . The left hand side will be referred to as the ACW and CW rotation average components for negative and positive sign respectively. It can be shown that the rotation average components in Eq. (14) are equivalent to the corresponding azimuthal wave components defined by the Bloch kernels ψ ±1 The proof of this is presented in Appendix A for completeness. The phase of the component in Eq. (15) is determined from the temporal phase ωt 0 selected for the first step in the rotation averaging process. Eq. (14) can also be used to obtain higher spatial harmonics of the HRR by changing the relation between the rotation term and the time step. In the current study the higher harmonics are of negligible order, and will not be considered [66] . The technique can also be slightly modified to work with different injector types in the same annular combustor, as long as the distribution of the different injectors is cyclic with more than a single period in the azimuthal direction.
Another way to consider this process is as follows. If the total phase average response corresponds to oscillations which travel in both CW and ACW directions around the annulus once per cycle, each different flame responds with a delay in phase to these oscillations. The rotational average essentially travels with the oscillations; averaging together the response of each individual flame at the same phase in the oscillation cycle. As there are both ACW and CW oscillations, two rotational averages can be calculated, yielding the average response of all flames to each phase in the oscillation cycle. In practice the phase average was divided into 3 N = 36 equally spaced phase bins, which can be thought of as 3 individual time series where the bins within each individual series is separated by t = 2 π / ( ωN ) . The start time t 0 of the 3 series are separated by a time delta of 2 π / ( 3 ωN ) , providing 3 times better interpolated temporal resolution compared to dividing it into N = 12 bins. For each of the 3 time series the rotation average in Eq. (14) is calculated for the real valued phase average images where the rotations are performed by physically rotating the coordinate system of the images for each time step.

Reconstructed heat release rate fluctuations
The rotation averaged components in Eq. (15) are the Bloch azimuthal wave components in each direction of the phase average, and are by definition describing the average flame response to a first azimuthal mode. This can be utilized to recreate the phase average in Eq. (12) without the flame to flame differences ( q a ) by performing the opposite rotation 1 q , rec a ( r, θ , t l ) = q a ( r, θ + 2 π l/N, t 0 ) where t l = t 0 + 2 π l/ ( ωN ) . The response of the reconstructed phase average effectively eliminates flame to flame differences from the phase average response for each forced state, enabling clearer interpretation of structures and dynamics. In practical terms the real part of the reconstruction is obtained as the superposition of the 1 Equation (12) without flame to flame differences q a is retrieved by inserting Eq. (15) into Eq. (16) . two real valued rotation averages rotated a fixed angle in opposite directions, where the angle ±2 π l/N corresponds to time step t l in the phase average. The same method can be used to reconstruct the time series of the individual Bloch wave components ψ ±1 e i ( ±θ + ωt ) by only including one of the right hand side terms in Eq. (16) without converting to frequency space.

Velocity fluctuations and the Azimuthal FDF
The same rotational averaging procedure and reconstruction can also be applied to the axial velocity perturbations evaluated at the dump plane. This assumes that the axial velocity perturbations can be expressed as Bloch waves corresponding to the first azimuthal mode. Since the induced acoustic mode in the injector tube is predominantly determined through pressure coupling with the azimuthal pressure mode [59] , the induced axial velocity perturbations, considered at the injector centre location, are determined by the local azimuthal pressure mode amplitude. For brevity the identical derivation will not be shown, but the results are equivalent. Since the axial velocity perturbations are calculated for the complex valued Fourier amplitudes by default, the rotation averaged component of the axial velocity amplitude at the dump plane can be obtained from ˆ u axial θ j Here the N * instrumented injectors are assumed to be equidistantly spaced in the azimuthal direction, meaning N * = N/ 2 = 6 in the current configuration. Eq. (17) describes the axial velocity perturbations as two counter propagating azimuthal waves, and the left hand side is therefore denoted an azimuthal axial velocity perturbation component. In a similar manner as Eq. (16) the reconstructed axial velocity Fourier amplitude can simply be obtained from ˆ u , rec axial θ j = ˆ u axial θ j As a shorthand, and to be more in line with the notation of the sector averages ( ·) j , the following notational convention is intro- One additional major advantage of the velocity signals as shown in Eq. (18) is that the axial velocity perturbations are obtained for all the injectors. This doubles the amount of usable data, as only every other injector is instrumented due to the cost of acquiring pressure signals and the required storage. It is then natural to define the Azimuthal Flame Describing Function as where the spatial averaging is again performed over one flame sector, as defined in Fig. 1 . The domain of the functions ˆ q ± and q sectors is a single flame sector, however due to imperfections in centering of the flame sector masks and small differences caused by defining the masks on a square grid, the Azimuthal FDF is calculated as the mean of all N flame sectors to reduce effects of this.

Azimuthal heat release rate mode
In the previous section the Azimuthal Flame Describing Function was introduced. This can be used to directly express the normalised, spatially sector integrated heat release rate fluctuations of each flame j as removing the small difference caused by the imperfect sector masks. The two rotation average components of the velocity have the following property according to the underlying Bloch theory where θ j = 2 π j/N. This means the expression in Eq. (21) describes two 1D counter propagating waves with analytical signal q , rec a j where the amplitudes Q ± are given by This formulation makes it straight forward to obtain the normalised, sector integrated HRR of a flame for a given pair of azimuthal velocity perturbations from a pair of Azimuthal FDFs that are either specific to the forced state or generalised for all the forced states. Using generalised FDFs enables easy interpolation between the different forced states, meaning an arbitrary pressure mode within the range of perturbation amplitudes of the Azimuthal FDFs can be modeled. Therefore, the pair of Azimuthal FDF components should be obtained for the range of amplitudes and mode orders of interest, either from experiments or high fidelity simulations. Eq. (23) means the HRR oscillations can be described as a HRR mode equivalent to the pressure mode in Eq. (5) . Therefore, the nature angle of the HRR mode is introduced as [18] As previously mentioned, the axial velocity perturbations are assumed to be distributed in the same manner as the pressure distribution in Eq. (5) , as the induced axial velocity is determined by the local pressure amplitude. Therefore, the nature angle of the pressure χ is the same as the nature angle calculated from the amplitudes of ˆ u axial ± . This implies the nature angle of the HRR mode χ q may not be the same as nature angle of the pressure mode χ in general. Assuming the Azimuthal FDFs are linear, the only case where χ q = χ for all χ is when the two Azimuthal FDFs have the same gain.

Results and discussion
In this section the conventional describing function ( F DF ) will be introduced in Section 3.1 to examine the response of the system through a metric similar to that applied in more conventional single flame configurations. This will also enable the study of the linearity of the system and quantify the spread in response of different flames. Then the difference in the response of the different flames are examined in Section 3.2 before studying the response at the pressure node location in Section 3.3 . The two Bloch kernels, obtained from the corresponding spinning forced states, are presented in Section 3.4 , showing the structure of the heat release rate oscillations. The new Azimuthal describing functions based on the Bloch wave response are presented in Section 3.5 . Finally the differences between the azimuthal pressure mode and the azimuthal HRR mode is studied detail in Section 3.6 . Figure 6 presents the heat release rate (HRR) fluctuations for all the forced states on the Poincaré sphere in Fig. 4 . Each forced state corresponds to 6 points in the figure, one for each instrumented injector, making a total of 738 points. The HRR fluctuations are the fluctuations of the spatially integrated heat release rate for each flame sector (shown in Fig. 1 b). The color of the points represents the nature angle χ of the forced states, and the spinning states are plotted in the foreground for improved viewing clarity. The points characterise the response of the flame through a similar measure as the conventional FDF. The main differences are the HRR fluctuations are plotted against velocity instead of frequency, allowing for evaluation of the linearity of the system, and the HRR is not normalised by the velocity, both similar to Balachandran et al. [23] . Plotting the F DF as separate HRR and velocity fluctuations enables the inclusion of points of relatively high HRR fluctuations for small velocity perturbations ˆ u axial , j / U bulk , such as the cluster denoted LPHR highlighted in Fig. 6 . The magnitude of the flame response to the anti-clockwise (ACW) spinning forced states is observed to generally be larger than the response to the clockwise (CW) spinning ones for all perturbation amplitudes, as demonstrated in [47] . The standing modes are mostly in between the outer limits created by the spinning forced states, with an interesting exception for the points in the Low Perturbation High Response ( LPHR ) cluster in Fig. 6 . This will be shown to depend on the location of the flame relative to the mode orientation θ 0 in Section 3.3 . Sections 3.5 and 3.6 will also give a plausible reason for why this difference in response exists. There is no clear difference in the trend of the ACW and CW forced states in the phase, which is observed to be mostly unaffected by the nature angle χ of the forced state. At low velocity perturbation levels, the spread in the phase is observed to increase, spanning the full range of phases at negligible perturbation levels. This is related to the LPHR cluster in Fig. 6 , and will be discussed more in relation to Fig. 9 .

Flame and pressure mode interaction
It is also instructive to normalise the HRR oscillations by the axial velocity oscillations in order to describe the gain parameter of the FDF. This is presented in Fig. 7 , and the system is observed to have a close to constant trend for the full range of axial velocity perturbations. The observation suggest the flame response can be considered linear, or close to linear, for the full range even at relatively high axial velocity perturbations. This is also consistent with the frequency being well into the gain cut-off range [11] . The observed approximate linearity of the system supports the use of the superposition assumption in Section 2.6 , which will be used in Section 3.4 and subsequent sections. Similarly to the observations in Fig. 6 , the ACW and CW spinning forced states are shown to have a different response to similar perturbation levels in Fig. 7 . However, the difference in the gain for the two spinning directions is increasing for decreasing perturbation amplitudes. Since the spread in the HRR response in Fig. 6 is similar for all perturbation amplitudes, the spread in the gain is inversely proportional to the perturbation amplitude in Fig. 7 .

Flame response of different injectors
The scatter observed in Figs. 6 and 7 is also influenced by small differences from flame to flame. In contrast to most other forced studies based on a single or a few flames, the current study employs N = 12 injectors which are identical by design, subjected to the same azimuthal pressure mode. The dominantly spinning forced states, here chosen to be | χ| ≥ π / 5 , therefore subject all the flames to approximately the same local pressure and axial velocity perturbations. This is not the case for the standing and mixed  Fig. 6 . Vertical axis is limited to aid readability, due to very large gains for some cases with very low axial perturbation levels. modes, where different flames experience different local pressure fluctuation amplitudes, resulting in different induced axial velocity perturbations. This makes the spinning forced states ideal to study the flame to flame differences, which is sho wn in Fig. 8 as a boxplot for the flames at the 6 instrumented injectors. The gain is calculated according to Eq. (7) based on the response for all forced states with the required nature angle, and the median value is marked by the solid horizontal lines within the box which it- self marks the range from the 25th to 75th percentile. Due to only modes of sufficient amplitude satisfying the nature angle requirement, the lowest axial velocity fluctuation amplitude is approximately ˆ u axial , j / U bulk ≈ 0 . 05 . Significant differences are observed between the different flames when subject to spinning waves travelling in both directions, suggesting there are flame to flame differences in the response to the same axial velocity perturbation signals for both forced state directions. In Section 3.5 these differences in response will be averaged out to ease interpretation of the average flame response. This is the first time that the level of loss of global rotational symmetry in a nominally symmetric annular combustor is assessed. It is then interesting to consider the work of [53] , where a very similar scatter in the values of the flames gain was considered (see their Fig. 5 ), for a self-excited, azimuthal thermoacoustic instability. They show that this level of flame to flame variation is responsible for a preferential orientation of the standing part of the acoustic pressure field as observed in experiments. They also show that this loss of global rotational symmetry does not however significantly affect the nature angle χ of the acoustic pressure field, i.e. whether the system spins or not [53, Fig. 6] . Therefore, it was concluded that for flame to flame differences of the order reported here, existing criteria for the existence and stability of standing and spinning self-excited solutions, obtained assuming rotational symmetry, still hold [67] .
It is also interesting to note the difference between the ACW and CW forced states is not the same for all the flames. Some flames show a significant difference, while others show a small difference in the median values. However, the mean gain for all the flames are higher for the ACW forced states compared to the CW forced states, further supporting the observation that there is a clear difference in the flame response to the direction of the forced states.

Response at the pressure node
The cluster of points from the standing states denoted LPHR in Fig. 6 is unexpected, as the points correspond to relatively high HRR response to very low axial velocity perturbations. To study the physical location of these flames relative the pressure mode, the data in Fig. 6 are plotted for the standing modes ( | χ| < π / 20 ) in Fig. 9 with the color determined by the distance from the pressure anti-node . The LPHR cluster is observed to consist of flames located approximately ±π / 2 away from the pressure mode antinode, which is the location of the pressure node. In the pressure node, the induced axial perturbations are negligible independent of the coupling impedance between the annular chamber and the injectors in the case of predominantly standing modes, | χ | < π/ 20 . Points are colored by the azimuthal distance from the pressure anti-node, and the gray points represent forced states | χ | ≥ π/ 20 for visual reference. The points with low perturbation and high HRR response (LPHR) correspond to the node of the pressure mode, where the induced axial velocity perturbations are small even for a high amplitude pressure mode. The corresponding phase is also observed to deviate from the expectations.
injector tube [59] , as observed in Fig. 9 . This is the likely cause of the greater spread in the phase, as the phase of the axial velocity in the FDF definition in Eq. (7) is less well defined at the azimuthal pressure node. However, for a given pressure mode the azimuthal velocity perturbations are at the maximum at the pressure node location for the standing modes [59] . Perfectly symmetric flames would not be affected by such velocity perturbations due to the symmetry of the flames [68] . When this symmetry is broken the transverse velocity perturbations have theoretically been shown to play a role in the response of the for a flame enclosed in a cylindrical can enclosure [50] . However, in Section 3.6 it will be argued that this effect is most likely a direct consequence of the systematic local symmetry breaking in the current configuration.

Bloch kernels
To study the difference in HRR structures in the response, the real components of the Bloch kernels ψ −1 and ψ +1 , obtained from an ACW and a CW forced state of similar amplitude, are shown for the zeroth flame in Fig. 10 . The real part was taken at the point where the mean phase of the spatially averaged Bloch kernel is zero, approximately corresponding to the structure at the point of maximum heat release rate. The most striking feature of the response is the large region of continuous high response, highlighted by the yellow dashed rectangle in the figure, changes from the right hand side of the flame for the ACW kernel ( top ) to the left hand side for the CW kernel ( bottom ). Another notable feature is that the propagation of heat release rate fluctuations from the center to the outer parts of the flame form a spiral-like pattern with regions of high and low response. This response can be contrasted with the response of the same injector when axially forced upstream in a single flame geometry, where the propagation of regions of high and low heat release rate was observed to form hexagonal rings propagating from the center of the flame [11] . While the change of side of the high heat release rate region can be attributed to the change in forcing direction, the spiral-like pattern of the propagation of the fluctuations means the response is not purely mirrored.
Adding the two kernels would result in a crescent-like shape around the flame along the outer wall and between adjacent injectors, which was observed for a self-excited standing mode in a similar setup [22] . The same study also concludes that the response is stronger on the side of the flame facing the incident acoustic wave for a 15 flame configuration, which is consistent with the observations made from the Bloch kernels in Fig. 10 . Dawson and Worth [22] speculated that this might be related to nonnegligible transverse velocity oscillations, which could cause vorticity cancellation similar to the process which occurs in a jet in crossflow during vortex formation [69] . The phenomenon was later observed in a single flame by Saurabh and Paschereit [51,52] , who also attributed it to a timing dependent interaction between the transverse and axial velocity perturbations. It is also worth noting the peak amplitude of regions of high HRR fluctuation amplitude in the figure is about 45 % higher for the ACW kernel ψ −1 than for the CW kernel ψ +1 . This might contribute to the lower response observed for the CW modes in Fig. 6 . Therefore, despite large differences in orientation, structural differences in the response appear to be slight. However, the response does show a significant reduction in heat release response amplitude, which is dependent on forcing direction. This phenomenon can also be observed directly from the phase average and rotation averages, which can be found in Appendix B in the supplementary material online.

Azimuthal flame describing functions
The conventional FDFs presented in Section 3.1 suggested the HRR response to the two different directions are different, but the picture is complicated by the fact that 6 flames are presented for each forced state. Each of these flames might be subjected to a  Fig. 6 . Both row presents the same data, with the component of interest presented in color and the data from the opposite component is shown in grey for visual reference. The solid black line is a linear regression through the origin of the magnitude of response, with the slope of the upper line (ACW) being 66 % larger than the lower one (CW). There are some differences in the phase at lower perturbation levels, but no significant difference is observed for the higher perturbation levels. different local pressure amplitude, and as a result different axial acoustic velocity fluctuations for a single forced state. Therefore, the Azimuthal Flame Describing Functions F DF ± , first introduced in Section 2.6 , are presented in Fig. 11 as a more suitable measure of the response of the system. Each forced state results in a single point for the ACW function ( F DF − , top row ) and a single point for the CW function ( F DF + , bottom row ). The ordering of the nature angles in the horizontal direction is observed to be inverted for the two cases. Considering the ACW component in the top row, the general trend is for the points to go from negative nature angles (CW mode) to positive nature angles (ACW mode). Similarly, for the CW component the nature angles go from positive on the left to negative on the right in general. This follows the expected behaviour due to the relative magnitude of the azimuthal wave components of the axial velocity perturbations being a function of nature angle χ.
Focusing on the magnitude of the response, there is a significant difference in the response of the two components. The ACW component, traveling with the bulk swirl in the outer half of the annulus, has a significantly ( 66 % ) higher response gain trend than the CW component for the same azimuthal axial velocity perturbation magnitude. The trend slope, representing the gain of the Azimuthal FDF, of both curves are approximately constant, with no significant scatter between the different forced states. The reduced scatter suggests the linearity assumption on the heat release rate in Section 2.6 holds. It should also make it possible to interpolate for all possible combinations of nature and orientation angles within the studied range of amplitudes through Eq. (23) . Another interesting observation is the lack of a cluster of points with high response for low perturbations, similar to LPHR in Fig. 6 . The cluster observed in the previous plots is in fact a direct effect of the difference in the Azimuthal FDF gain observed in Fig. 11 , but this will be discussed in detail in Section 3.6 . The reconstruction of the HRR and velocity fluctuations from this linear decomposition are also shown to represent the conventional F DF s well in Appendix C found in the supplementary material online, with the same cluster of points retained.
The phase of the two azimuthal components, presented in the right column of Fig. 11 , exhibit a very similar trend at the highest axial velocity perturbations. At smaller perturbation levels there are some differences in the phase between the two components, suggesting there are small timing differences in the response. However, both components show a relatively small spread compared to the phase observed for standing modes in Fig. 6 , which will be discussed in more detail in Section 3.6 . For negligible perturbation levels the strongly spinning modes have a few points with random scatter in the phase, but this is most likely an artefact of the very low perturbation amplitude which makes the phase more susceptible to noise.
The above observations of Fig. 11 show a clear case of the first azimuthal HRR mode splitting into two distinct modes. For an annular system with a mean swirl it has recently been shown that the first azimuthal acoustic mode will split in both growth rate and frequency [43,45] . Due to the low flow rate, and large cross section of the combustion chamber in the current study, the induced net swirl gives rise to a flow of low Mach number [70] . From simple geometric considerations the azimuthal bulk velocity in the inner and outer part of the annulus induced by the swirl is estimated to be of the order of 1 m / s . Together with the forcing, this means the frequency splitting should be of negligible order. However, the HRR shows a significant splitting of gain of the two azimuthal directions. While a difference in frequency of the acoustic response has been observed experimentally previously [44] , the effect on the HRR mode has not been observed systematically before due to a lack of systematic studies of forced annular combustors. The observation that the Azimuthal F DF s have a significant difference in gain has several implications for the response, some of which will be studied in Section 3.6 .

Effect of local symmetry breaking
The difference in response for the ACW and CW directions is most likely caused by the local symmetry breaking introduced by having swirling flames in an annular combustor of finite width and radius. Considering a perfect geometry, the combustor has reflectional symmetry when the flames are without swirl, making it impossible to distinguish a mirrored system from the original at any flame location, as observed for a practical geometry in [21] . Therefore, the response has to be the same for both spinning directions from a symmetry point of view. Similarly, for setups with a local 180 • rotational symmetry around the flame center, such as conventional transverse forcing setups [35][36][37][38][39] , it is also impossible to separate the two forcing directions, suggesting the response in both directions has to be the same. Since the argument relies on a rotational symmetry, it is independent of the swirl number S of each injector. In an annular combustor of finite width and radius with swirling flames both of these symmetries are broken, allowing for a splitting of the HRR mode. The local reflectional symmetry is broken by the swirl, and the local rotational symmetry around the centre of a flame is broken by the finite width and curvature of the combustion chamber.
In practice, small inherent asymmetries in single flames have been observed to break the theoretical equivalence of the two forcing directions in a novel combined axial and transverse forcing setup [51,52] . This is supported by the analytical considerations of Acharya and Lieuwen [50] , which showed the response can become sensitive to the transverse velocity fluctuations when the flame is not fully symmetric. However, in the current study there are N = 12 different flames, which are combined together to produce an "average flame" through rotation averaging. Assuming the inherent asymmetries in a single flame are due to small manufacturing and mounting imperfections, the small asymmetries should not be systematically distributed, which is supported by the difference from flame to flame observed in Fig. 5 . This suggests that the impact of the individual flame asymmetries is reduced compared to a single flame setup, and it is conjectured that the main symmetry breaking effects stem from the systematic breaking of the local symmetry by the geometry combined with swirling flames.
To quantify the degree of asymmetry introduced by the swirl and finite dimensions of the combustion chamber, two nondimensional numbers are proposed. The first number is the conventional swirl number S of each injector. The effect of the geometry is proposed to be quantified through a new non-dimensional curvature number The non-dimensional curvature number , which is bounded in [0,1], describes the width of the annular combustion chamber relative to the mean radius. The qualitative relation between the transferability of the conventional FDF and the two asymmetry parameters S and is presented in Fig. 12 . The conventional FDF is expected to be highly transferable for combustors without swirl, or for geometries with = 0 , corresponding to either infinitely thin geometries or a finite width annulus of infinite radius, which is equivalent to conventional transverse forcing setups [35][36][37][38][39] . To experience a difference in the magnitude of the Azimuthal FDF components, such as the difference observed in Fig. 11 , both S and have to be non-zero ( S ≈ 0 . 65 and = 0 . 24 in this case).

Pressure and HRR mode comparison
One important implication of the difference in gain between the two azimuthal flame describing functions observed in Fig. 11 can be illustrated by the sector integrated heat release rate time series compared to the pressure time series. Fig. 13 shows the HRR fluctuations for one period against the spatial location in the left column, in addition to the corresponding pressure mode shown both as a function of time for select locations ( top right ) and as a function of azimuthal angle θ ( bottom right ). The top left shows the normalised HRR directly from the phase average, while the bottom left shows the same based on the reconstructed version of the phase average. For curves based on the reconstruction, the envelope has the shape of a sinusoidal standing mode, with a pair of nodes and a pair of anti-nodes. The phase average has the same For injectors with no swirl and for very thin annuli the effect is absent: the flame response is the same in an annular combustor and in a single-flame test setup that allows transverse forcing in the orthogonal direction (assuming no flame-flame interaction occurs). When both these conditions are violated, the flame response to anti-clockwise and clockwise acoustic waves may differ, as found experimentally in this paper, and transferability is compromised. This is the case for most annular combustors, which have 0 < < 1 and swirling flames. general trend, but it is a bit harder to identify as a standing mode due to the flame to flame differences. Some locations have a response either higher or lower than what is expected of a sinusoidal standing mode, but fitting Eq. (23) results in the same HRR nature angle χ q ≈ 0 suggesting a standing mode. Therefore, the average dynamics of the system are easier to interpret based on the reconstructed average.
The most interesting point about Fig. 13 however comes from comparing the HRR mode, on the left, to the pressure mode shown on the right. The HRR mode was observed to be a standing mode, but the prescribed pressure mode has a relatively strong CW spin- The fact that the HRR nature angle χ q is significantly higher than the pressure mode nature angle χ is not only true for this single forced state. The HRR nature angle χ q and the pressure nature angle χ for all the forced states are shown in Fig. 14 , showing the nature angle of the HRR and pressure modes are different in general.
The nature angle of the HRR mode is typically higher than the nature angle of the pressure mode, with a decreasing difference at the physical limits of the pressure nature angle. Assuming the response can be approximated by the superposition of the two Azimuthal FDFs, as described in Eq. (23) , this can be predicted from the different slopes in Fig. 11 by inserting into Eq. (25) . The prediction based on the fitted gain in Fig. 11    angles χq and χ are the same. The solid black line is the predicted value of χq given a pressure mode nature angle χ , based on the linear behaviour described in Section 2.7 and the 2 linear trends of Fig. 11 . The observed χq follows the predicted line very well, with the largest scatter away from the line observed for very small pressure amplitudes which are more susceptible to noise. direction compared to the ACW direction ( | F DF − | < | F DF + | ) , the relation would be reversed, with χ q ≤ χ.
The points in Fig. 14 , which are based on the conventional phase average (with no azimuthal decomposition by rotation averaging), follow the predicted line well. There is some scatter, but it is mostly for points of lower pressure amplitude, which are more susceptible to noise. The slope of the predicted line is steeper than unity close to χ = −π / 4 , making the HRR nature angle χ q more sensitive to changes in the nature angle of the pressure mode χ. In the other end, close to χ = π / 4 , the slope of the predicted line is below unity, suggesting the HRR nature angle is less sensitive to changes in the nature angle of the pressure mode. Since the HRR nature angle χ q is predicted to generally be higher than χ and there is a hard limit of χ q ≤ π / 4 , the difference in nature angles cannot be large due to noise or other experimental inaccuracies close to χ = π / 4 . This is not the case close to χ = −π / 4 , and therefore larger scatter can be expected and is observed. The notion that the nature angle of the pressure mode and the HRR mode are not the same, or very similar, is somewhat unexpected and very interesting. This has never been observed before, as previous forcing studies have not had enough data or the right geometry to observe this, and the effect is not captured in models using the conventional FDF obtained in axially perturbed single flame setups. Interestingly, Fig. 11 b in [14] show some evidence of the nature angle of the HRR being higher than the nature angle of the pressure mode, although this was not mentioned at the time as the self-excited study lacked the necessary systematic exploration of nature angles required to confirm this behaviour. This is for a slightly different combustor, equipped with 18 injectors, compared to the 12 used in the current study.
In addition to studying the nature angle of the HRR and pressure mode, it is interesting to study the orientation angle of the standing component of the two modes. This is presented in Fig. 15 for the same forced states as Fig. 14 . For sufficiently high pressure amplitudes, and predominantly standing modes where the orientation angle is well defined, the orientation is observed to be the same for both the HRR and pressure mode. A small difference is observed when the amplitude is decreased slightly, but most points still stay within the azimuthal width of one flame sector. This might suggest that the HRR mode does not follow the pressure mode orientation as well for lower amplitudes. Some points of low amplitude, and the points close to the extreme nature angles χ = ±π / 4 show a relatively large difference in the orientation of the mode. However, in most of these cases the amplitude A st of the standing part of the pressure mode is below 100 Pa , with some points up to A st ≈ 275 Pa . This means all these points have a relatively weak standing component, and may therefore be more susceptible to noise and are less likely to perfectly follow to the pressure mode orientation.
The higher nature angle of the HRR mode compared to the pressure mode is most likely the main contributing factor to the non-negligible HRR observed for the standing modes in the LPHR cluster in Figs. 6 and 9 . This is exemplified by the HRR ( top ) and pressure ( bottom ) modes shown in Fig. 16 (which are the same as the forced standing state in Fig. B.1 (middle row) in the supplementary material online). Here, the nature angle of the pressure mode is χ = −0 . 01 π / 4 , while the nature angle of the HRR mode is χ q = 0 . 35 π / 4 . The axial velocity perturbations are determined by the local pressure mode, and are therefore negligible in the pressure node [59] shown in the bottom row of Fig. 16 . However, due to the significantly higher nature angle of the HRR mode, a finite amplitude response at the location of the pressure node is observed in the top plot of Fig. 16 . This is causing the relatively large HRR response to very low axial velocity perturbations, as observed in the LPHR cluster in Fig. 6 . In that figure the phase is also observed to have a much larger scatter for these low perturbation high response points. However, since the axial velocity is of negligible order, the phase is very sensitive to noise, as it is less well defined, and the scatter should be considered to be an artefact of the low perturbation amplitude in the specific injector. Generally, as observed in Fig. 9 , the points close to the node of a standing mode would be expected to have a higher response than that predicted by a conventional single flame setup due to the mixed mode nature of the corresponding HRR mode.
Therefore, for modeling it would make sense to use the azimuthal wave response, as defined in Eq. (23) from the Bloch theory approach. This formalism captures the response at the different locations in the annular geometry and enables the interpolation to an arbitrary pressure mode, up to the perturbation limit set by the physical limit of the forcing array. It also captures the effects of the finite curvature and width of the annular combustion chamber with swirling flames, like the difference in nature angle between the HRR and pressure modes, allowing it to be used in equations that assume infinitely thin combustion chambers. This is much simpler than using the full, conventional describing function for each different flame and every forced state. Therefore, it is suggested to use the Azimuthal FDF over the conventional FDF in models where the full frequency depdence is not required, as the conventional FDF result can be obtained by imposing an equal magnitude on each of the two Azimuthal FDF components.

Conclusions
The current study is the first to present the experimentally determined heat release rate response of an annular combustor with swirl for the full range of state space parameters of the first azimuthal mode up to a finite amplitude limit. This was achieved using azimuthal acoustic forcing, resulting in a 40 fold increase in the number of studied azimuthal forced states compared to previous studies in similar configurations. The response was first studied in terms of the conventional Flame Describing Functions (FDFs), based on acoustic axial velocity at the dump plane as input and heat release rate fluctuations as output, which suggested there is a difference between the different directions of the dominant first azimuthal component of the pressure mode. Additionally, it was observed that flames placed exactly at the pressure node location exhibit a relatively large heat release rate response despite negligible axial velocity perturbations in the case of standing forced states. This work has also assessed for the first time the effect that small flame to flame differences have on the respective flame transfer functions. This is the first time that this departure from the nominally rotational symmetric configuration, which is the reference in the design process, testing and modelling, is quantitatively assessed. For the annular combustor considered, the standard deviation of the flames' median gains is approximately 0.1 of their mean value for each spinning direction. Up to this value, we discuss how this level of global rotational symmetry does not impact significantly the nature angle χ of the acoustic field. In the industrial setting this number proves as a useful reference, for example when accepted tolerances on flame to flame differences are considered.
To reduce the complexity in the response the flame to flame differences were removed by the use of rotation averaging. The rotational averaging was shown to correspond to extracting the first azimuthal mode components of the response following Bloch theory. Bloch theory is used to remove flame to flame differences by defining the mean flame response from the known heat release rate of all the flames, subject to a forced state that is azimuthal. This also made it possible to reduce the full time series into two Bloch kernels, which determine the full, spatial response of the configuration in the absence of flame to flame differences to an azimuthal mode. The structure of two kernels was observed to be similar, but the peak amplitude of the anti-clockwise kernel was observed to be significantly higher (approximately 45 % ).
The complicated picture described by the conventional FDFs was simplified by introducing the concept of Azimuthal FDFs ( F DF ± ), based on the decomposition into azimuthal components by rotation averaging. In this framework the response of the two different spinning components was shown to be approximately linear for the full range of velocity amplitudes in the current study. A significantly higher gain ( 66 % ) was observed for the anticlockwise component. This is a direct result of the systematic symmetry breaking associated with an annular geometry combined with swirling flames, as in this formalism flame to flame differences are averaged out.
The difference in response of the two Azimuthal FDF components was shown to have several implications on the response, with the most prominent effect being the nature angle of the heat release rate mode χ q is in general different from the nature angle of the pressure mode χ = χ q in the current configuration. This again was shown to result in relatively large heat release rate response despite negligible acoustic axial velocity perturbations at the pressure node in the case of a perfectly standing pressure mode. Therefore, the effect of nature angle on the flame describing function should be taken into account when modelling annular combustion chambers with swirl. All the effects are described well by the Azimuthal FDF components, suggesting the characterization of the response of annular combustors to the first azimuthal mode is best described by this formalism compared to the more conventional FDFs. By systematically studying the response for the full range of state space parameters it was observed that the Azimuthal FDF components are not dependent on the nature angle of the pressure mode, making it possible to use this to interpolate between all possible state space combinations up to an amplitude limit of validity.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. ψ ±1 ( r, θ ∓ 2 π l/N ) e ±i ( θ ∓2 π l/N ) e i ω ( t+2 π l/ ( ωN ) ) , (A.4) where the periodicity is exploited in the last step by writing ψ ±1 ( r, θ ) = ψ ±1 ( r, θ ∓ 2 π l/N ) . This expression have the same general form as the rotation averaging definition in Eq. (14) . The last step is to insert Eqs. (12) and (A.3) into the rotation averaging definition in Eq. (14) to obtain the final expression q a ( r, θ , t 0 ) ± = ψ ±1 ( r, θ ) e ±i θ e i ωt 0 . It is worth noting all the steps in the rotation averaging process are linear operations. The only operations that are required are rotation, addition and division by a constant factor N. Therefore, the real and imaginary part of the analytical expression for the HRR q a can be treated separately. In practice the real part of the rotation average components is therefore calculated directly from the real valued phase averaged HRR images.

Supplementary material
Supplementary material associated with this article can be found, in the online version, at doi: 10.1016/j.combustflame.2021. 111565 .