Davide Laera
email: davide.laera@poliba.it

Thierry Schuller

Kevin Prieur

Daniel Durox

Sergio M Camporeale

Sébastien Candel

Flame Describing Function analysis of spinning and standing modes in an annular combustor and comparison with experiments

Keywords:

This article reports a numerical analysis of combustion instabilities coupled by a spinning mode or a standing mode in an annular combustor. The method combines an iterative algorithm involving a Helmholtz solver with the Flame Describing Function (FDF) framework. This is applied to azimuthal acoustic coupling with combustion dynamics and is used to perform a weakly nonlinear stability analysis yielding the system response trajectory in the frequency-growth rate plane until a limit cycle condition is reached. Two scenarios for mode type selection are tentatively proposed. The first is based on an analysis of the frequency growth rate trajectories of the system for different initial solutions. The second considers the stability of the solutions at limit cycle. It is concluded that a criterion combining the stability analysis at the limit cycle with the trajectory analysis might best define the mode type at the limit cycle. Simulations are compared with experiments carried out on the MICCA test facility equipped with 16 matrix burners. Each burner response is represented by means of a global FDF and it is considered that the spacing between burners is such that coupling with the mode takes place without mutual interactions between adjacent burning regions. Depending on the nature of the mode being considered, two hypotheses are made for the FDFs of the burners. When instabilities are coupled by a spinning mode, each burner features the same velocity fluctuation level implying that the complex FDF values are the same for all burners. In case of a standing mode, the sixteen burners feature different velocity fluctuation amplitudes depending on their relative position with respect to the pressure nodal line. Simulations retrieve the spinning or standing nature of the self-sustained mode that were identified in the experiments both in the plenum and in the combustion chamber. The frequency and amplitude of velocity fluctuations predicted at limit cycle are used to reconstruct time resolved pressure fluctuations in the plenum and chamber and heat release rate fluctuations at two locations. For the pressure fluctuations, the analysis provides a suitable estimate of the limit cycle oscillation and suitably retrieves experimental data recorded in the MICCA setup and in particular reflects the difference in amplitude levels observed in these two cavities. Differences in measured and predicted amplitudes appear for the heat release rate fluctuations. Their amplitude is found to be directly linked to the rapid change in the FDF gain as the velocity fluctuation level reaches large amplitudes corresponding to the limit cycle, underlying the need of FDF information at high modulation amplitudes.

Introduction

Many recent studies focus on combustion instabilities coupled by azimuthal modes in annular systems. There are yet few comparisons between predictions and well controlled experiments. The present investigation aims at filling this gap by developing a nu-• Two different operating conditions are considered, one leading to a spinning limit cycle and another one leading to a standing limit cycle. • This framework is then used to calculate the limit cycles of standing and spinning solutions and to compare the calculated oscillation with measurements on a laboratory scale test facility "MICCA" developed at EM2C laboratory. The amplitudes and phase relationships of pressure in the plenum and chamber and the heat release rate signals are compared for two different operating conditions. • Finally, two scenarios are tentatively proposed to explain the mode type selection. The first considers that the frequency and growth rate trajectories of initially spinning and standing modes determine the solution at the limit cycle. The second suggested by a reviewer considers that it is the stability of the limit cycle which determines the observed oscillation.

At this point, it is worth briefly reviewing some recent investigations of instabilities in annular devices. Combustion instabilities coupled by azimuthal modes are often studied by theoretical or numerical means. On the numerical level one finds a growing number of massively parallel large eddy simulations of annular chambers with multiple burners [START_REF] Staffelbach | Large eddy simulation of self excited azimuthal modes in annular combustors[END_REF][START_REF] Fureby | Les of a multi-burner annular gas turbine combustor[END_REF][START_REF] Wolf | Using les to study reacting flows and instabilities in annular combustion chambers[END_REF] . These calculations retrieve features of azimuthally coupled combustion instabilities observed in experiments in engine-like conditions [START_REF] Worth | Self-excited circumferential instabilities in a model annular gas turbine combustor: global flame dynamics[END_REF] . However, the complexity of the situation precludes direct comparison between calculations and observations.

One difficulty in performing such comparisons lies in the definition of the flame model. The first investigations of combustion instabilities coupled to azimuthal modes were performed by combining a low-order acoustic model of an annular combustor with a time delayed n -τ flame response (see for example [START_REF] Stow | Thermoacoustic oscillations in an annular combustor[END_REF][START_REF] Evesque | Low-order acoustic modelling for annular combustors: validation and inclusion of modal coupling[END_REF] ). This kind of model is also assumed in the recent analytical studies developed in [START_REF] Parmentier | A simple analytical model to study and control azimuthal instabilities in annular combustion chambers[END_REF][START_REF] Bauerheim | An analytical model for azimuthal thermoacoustic modes in an annular chamber fed by an annular plenum[END_REF] to analyze the linear stability of annular chambers fed by an annular plenum with multiple discrete burners. Both spinning and standing modes are predicted depending on the circumferential symmetry of the system. Circumferential instabilities of industrial combustors were analyzed in [START_REF] Krebs | Thermoacoustic stability chart for high-intensity gas turbine combustion systems[END_REF][START_REF] Campa | Prediction of the thermoacoustic combustion instabilities in practical annular combustors[END_REF] by means of a Helmholtz solver approach. In these studies, the flame response was modeled by a n -τ description with parameters retrieved from CFD calculations of the steady combustion process. These previous studies carried out with linear tools could not account for finite amplitude effects that determine the oscillation frequency and level at the limit cycle. These features were considered for example in [START_REF] Pankiewitz | Time domain simulation of combustion instabilities in annular combustors[END_REF][START_REF] Stow | A time-domain network model for nonlinear thermoacoustic oscillations[END_REF] by combining different numerical strategies for acoustic propagation with a nonlinear flame model in the time domain. Numerical control strategies for annular configurations featuring spinning limit cycles were developed in [START_REF] Morgans | Model-based control of combustion instabilities in annular combustors[END_REF][START_REF] Illingworth | Adaptive feedback control of combustion instability in annular combustors[END_REF] . Spinning and standing modes were observed in [START_REF] Pankiewitz | Time domain simulation of combustion instabilities in annular combustors[END_REF] depending on whether the mean flow velocity was neglected or considered in the simulations. For a circumferential instability in an axisymmetric geometry, the spinning waveform was always preferred to the standing mode in the simulations carried out in [START_REF] Stow | Thermoacoustic oscillations in an annular combustor[END_REF] . However, no comparison with experiments were reported in these works.

On a theoretical level, criteria for appearance of spinning or standing modes have been derived by considering the dynamics of azimuthal modes coupled by a nonlinear flame model expressed in terms of pressure perturbations alone. One of the first analysis was carried out in [START_REF] Schuermans | Non-linear combustion instabilities in annular gas-turbine combustors[END_REF] with a saturation function linking heat release and pressure fluctuations. In this framework, the dynamics of azimuthal modes is described by two harmonic oscillators which are nonlinearly coupled. The stability of standing and traveling waves at limit cycle can then be assessed. A further simplification was later introduced by assuming that the system behaves like a Van der Pol oscillator [START_REF] Noiray | Investigation of azimuthal staging concepts in annular gas turbines[END_REF][START_REF] Ghirardo | Azimuthal instabilities in annular combustors: standing and spinning modes[END_REF][START_REF] Noiray | On the dynamic nature of azimuthal thermoacoustic modes in annular gas turbine combustion chambers[END_REF] . Results indicate that the spinning or standing nature of the unstable mode originates from the nonlinearity and non-uniformity of the flame response and can be influenced by different factors, like transverse velocity fluctuations [START_REF] Dawson | Flame dynamics and unsteady heat release rate of self-excited azimuthal modes in an annular combustor[END_REF][START_REF] Oconnor | Recirculation zone dynamics of a transversely excited swirl flow and flame[END_REF] or turbulence, which can stochastically disturb the limit cycle amplitudes. This nonlinear flame model was used for example in [START_REF] Bothien | Analysis of azimuthal thermo-acoustic modes in annular gas turbine combustion chambers[END_REF] to reproduce the dynamical behavior observed in a real engine. Another comparison is presented in [START_REF] Noiray | Deterministic quantities characterizing noise driven hopf bifurcations in gas turbine combustors[END_REF] between numerical and experimental growth rates calculated by means of a system identification technique, but the oscillation levels of the different pressure signals are not shown.

Recently, Ghirardo et al. [START_REF] Ghirardo | State-space realization of a describing function[END_REF][START_REF] Ghirardo | Weakly nonlinear analysis of thermoacoustic instabilities in annular combustors[END_REF] managed to introduce in their time domain model a more reliable FDF, linking heat release and pressure disturbances by a time-invariant nonlinear operator. Criteria for self-sustained thermo-acoustic instabilities coupled by spinning and standing modes were then derived by examining the stability of the analytical solutions at limit cycles. This framework was then tested with experimental data from Bourgouin et al. [START_REF] Bourgouin | Characterization and modeling of a spinning thermoacoustic instability in an annular combustor equipped with multiple matrix injectors[END_REF] where a stable spinning mode is observed at limit cycle. Their analysis confirmed that for the operating condition explored, there was a stable spinning solution and that standing solutions, if they existed, were unstable. In all of these previous studies there are no direct comparisons between predictions and measurements for pressure and heat release rate oscillations and only limited validations of model predictions for different operating conditions.

The difficulties encountered in these various investigations are compounded by the presence of multiple flames which respond collectively over a wide frequency range, and by the modal density in the annular geometry when the size of the system is large like in industrial gas turbine combustors. One possible simplification consists in considering that the heat release from the different burners is uniformly distributed over the circumference of the annular chamber. Following this approach, Bourgouin et al. [START_REF] Bourgouin | Characterization and modeling of a spinning thermoacoustic instability in an annular combustor equipped with multiple matrix injectors[END_REF] developed an analytical one-dimensional framework to represent the dynamics of the laboratory scale MICCA annular combustor. Assuming a simplified flame response, the spinning instability recorded during experiments was reproduced in terms of frequency and amplitude of velocity fluctuations at the limit cycle. A theoretical interpretation was also proposed for the angular shift observed in the experiments between the nodal lines in the plenum and in the chamber. However, this analysis was carried out for a fixed frequency and fixed oscillation level corresponding to the values observed in the experiment at the limit cycle.

On the experimental level there are relatively few data sets corresponding to instrumented conditions that can be used to benchmark models and simulations. Most of the measurements performed on real systems consist of unsteady pressure signals with no access to the flame dynamics [START_REF] Krebs | Thermoacoustic stability chart for high-intensity gas turbine combustion systems[END_REF][START_REF] Bothien | Analysis of azimuthal thermo-acoustic modes in annular gas turbine combustion chambers[END_REF][START_REF] Mastrovito | Analysis of pressure oscillations data in gas turbine annular combustion chamber equipped with passive damper[END_REF] . Both spinning and standing mode patterns were observed in the laboratory scale annular device equipped with low swirl injectors that was developed in the engineering department of Cambridge University [START_REF] Worth | Self-excited circumferential instabilities in a model annular gas turbine combustor: global flame dynamics[END_REF][START_REF] Worth | Modal dynamics of self-excited azimuthal instabilities in an annular combustion chamber[END_REF] . This setup allows heat release rate measurements and flame dynamics analysis through optical windows, but the flame transfer function was not determined and pressure signals were not recorded in the combustion chamber.

This article is organized as follows. A novel procedure is derived in Section 2 that combines a Helmholtz solver with sixteen independent FDFs. This is used to determine the limit cycle conditions by means of a weakly nonlinear stability analysis. Depending on the nature of the mode being considered, the FDFs are assumed to operate at equal or at different velocity fluctuation levels for the examination of the dynamics of spinning and standing modes, respectively. This numerical procedure is validated in Appendix A in an idealized geometrical configuration with a simplified flame model by retrieving the amplitude and stability properties of the theoretical limit cycles [START_REF] Noiray | Investigation of azimuthal staging concepts in annular gas turbines[END_REF] . Section 3 describes the MICCA combustor experiment with 16 matrix laminar injectors and the numerical framework used for the stability analysis. The FDF determined for one of the matrix injectors is presented in Section 4 for the first operating condition leading to a spinning mode limit cycle in the MICCA combustor. The same type of analysis is repeated in Section 5 for the second operating point leading to a standing mode limit cycle. The frequency and amplitude of velocity fluctuations predicted at limit cycle are used to reconstruct the pressure oscillations in the plenum and in the chamber and the heat release rate fluctuation signal. These signals are then compared with microphone records and photomultiplier measurements at two different flame locations. A mode type selection analysis is presented in Section 6 . For the two operating points, frequencygrowth rate trajectories are calculated and the stability properties of the simulated limit cycles are discussed. These analyses are tentatively used to determine the limit cycle structure.

Weakly nonlinear stability numerical analysis

The methodology used to assess the stability of the annular combustor acoustic modes and determine the limit cycle state of the system follows that developed for the FDF based weakly nonlinear stability analysis of single burner setups [START_REF] Noiray | A unified framework for nonlinear combustion instability analysis based on the flame describing function[END_REF][START_REF] Palies | Nonlinear combustion instability analysis based on the flame describing function applied to turbulent premixed swirling flames[END_REF][START_REF] Silva | Combining a Helmholtz solver with the flame describing function to assess combustion instability in a premixed swirled combustor[END_REF][START_REF] Han | Prediction of combustion instability limit cycle oscillations by combining flame describing function simulations with a thermoacoustic network model[END_REF] . In the time domain one has to consider a wave equation including a damping term defined by a first order time derivative of the pressure fluctuations p multiplied by a damping rate δ ( s -1 ) [START_REF] Kinsler | Fundamentals of acoustics[END_REF] and a source term representing effects of the unsteady heat release rate disturbances ˙ q :

∂ 2 p ∂t 2 + 4 πδ ∂ p ∂t -∇ • c 2 ∇ p = (γ -1) ∂ ˙ q ∂t , (1)
where it has been assumed that the mean pressure p is essentially constant so that ρc 2 is also constant. Assuming that all fluctuations are harmonic x = ˆ

x exp (-iωt) , one obtains the following equation in the frequency domain:

ω 2 c 2 ˆ p + i ω4 πδ c 2 ˆ p + ρ∇ • 1 ρ ∇ ˆ p = i γ -1 c 2 ω ˆ ˙ q, ( 2
)
where ω denotes the complex angular frequency. The mean density ρ, speed of sound c and specific heat ratio γ distributions are specified. In this frequency domain, the analysis is carried out by using the finite element method (FEM) based on the commercial software COMSOL Multiphysics . This code solves the classical Helmholtz equation in which heat release rate fluctuations ˆ ˙ q are treated as pressure sources.

A nonlinear description of the interaction between combustion and acoustics is required to capture the limit cycles of a thermoacoustic system. If the flame is compact with respect to the wavelength of the unstable mode, the dynamics of the flame may be represented in terms of a global FDF F [START_REF] Silva | Combining a Helmholtz solver with the flame describing function to assess combustion instability in a premixed swirled combustor[END_REF] , where the FDF gain and phase lag are function of the amplitude of the incoming perturbation [START_REF] Noiray | A unified framework for nonlinear combustion instability analysis based on the flame describing function[END_REF][START_REF] Han | Simulation of the flame describing function of a turbulent premixed flame using an open-source les solver[END_REF] . In the frequency domain F links relative heat

release rate fluctuations ˆ ˙ Q/ ˙ Q to relative velocity fluctuations | ˆ u / u |
measured at a reference point of the system. The FDF is a complex function expressed in terms of a gain G and phase ϕ as follows:

F (ω r , | u / u | ) = ˆ ˙ Q (ω r , | u / u | ) / Q | ˆ u / u | = G (ω r , | u / u | ) exp iϕ(ω r , | u / u | ) , (3)
where | u / u | stands for the relative velocity fluctuation level, with u the root-mean-square of the velocity signals taken at the reference position in the injector unit j and ω r corresponds to the real part of the complex frequency ω. A weakly nonlinear approach is used to couple Eq. ( 3) with Eq. ( 2) retrieving the solution of the nonlinear problem as a perturbation of a linear problem. This is achieved by linearizing the FDF by fixing a velocity fluctuation level | u / u | . The finite element discretization of this set of linearized equations along with the boundary conditions results in the following eigenvalue problem [START_REF] Nicoud | Acoustic modes in combustors with complex impedances and multidimensional active flames[END_REF][START_REF] Laera | A finite element method for a weakly nonlinear dynamic analysis and bifurcation tracking of thermo-acoustic instability in longitudinal and annular combustors[END_REF] :

[ A ] P + ω [ B (ω) ] P + ω 2 [ C ] P = [ D (ω) ] P , ( 4
)
where P is the pressure eigenmodes vector. The matrices [ A ] and [ C ] contain coefficients originating from the discretization of the Helmholtz equation, [ B ( ω)] is the matrix of the boundary conditions and of the damping and [ D ( ω)] represents the source term due to the unsteady heat release. With the introduction of the heat release rate the eigenvalue problem defined by Eq. ( 4) becomes nonlinear and is solved with an iterative algorithm. At the k th iteration Eq. ( 4) is first reduced to a linear eigenvalue problem around a specific frequency k :

( [ A ] + k [ B ( k ) ] -[ D ( k ) ] ) P + ω 2 k [ C ] P = 0 , ( 5
)
where k = ω k -1 is the previous iteration result. The software uses the ARPACK numerical routine for large-scale eigenvalue problems. This is based on a variant of the Arnoldi algorithm, called the implicit restarted Arnoldi method [START_REF] Lehoucq | ARPACK users' guide: solution of large-scale eigenvalue problems with implicitly restarted Arnoldi methods[END_REF] . This procedure is iterated until the error defined by = | ω kk | becomes lower than a specific value, typically 10 -6 . Once convergence is achieved, the real part of ω yields the oscillation frequency, f = -( ω)/2 π Hz, while the imaginary part of ω corresponds to the growth rate α = (ω) / 2 π s -1 that allows the identification of unstable modes:

p ∝ exp ( 2 παt -i 2 π f t ) .
If α is positive, the acoustic mode is unstable and the amplitude of fluctuations grows with time. If α is negative, the acoustic mode is stable and perturbations decay with time. The eigenvalue procedure is repeated by incrementing the amplitude level until a limit cycle condition is reached when

α = 0 .
Differently from Silva et al. [START_REF] Silva | Combining a Helmholtz solver with the flame describing function to assess combustion instability in a premixed swirled combustor[END_REF] approach for single burner setups, in a multi-burner annular combustor the linearization is performed for each of the FDFs assumed in the model. In the present numerical framework, the distribution of velocity fluctuation levels between the FDFs is fixed a priori depending on the spinning or standing nature of the azimuthal mode under investigation:

• In a spinning mode, the nodal line rotates at the speed of sound, however, the oscillation amplitude is uniform all around the chamber and the velocity fluctuation level | u / u | is the same for each burner. In this case, the FDFs of the individual burners have the same complex value [START_REF] Bourgouin | Characterization and modeling of a spinning thermoacoustic instability in an annular combustor equipped with multiple matrix injectors[END_REF] . Mathematically this is formulated as follows:

| u / u | ( ) = C. ( 6
)
In this expression, is a vector containing the angular coordinates of the reference points where the velocity fluctuation is specified allowing the calculation of the FDF. Each reference point lies on the injector axis 20 mm below the injector outlet.

• In a standing mode, each injector operates with a different amplitude of oscillation depending on its relative position with respect to the nodal line [START_REF] Durox | Nonlinear interactions in combustion instabilities coupled by azimuthal acoustic modes[END_REF] . As a consequence, different am-

| u / u | ( ) = ψ ( ) /ψ max | u / u | j , (7)
where ψ ( )/ ψ max is the normalized azimuthal eigenmode structure. The numerical implementation of this model in the Helmholtz solver framework is not trivial because the mode ψ is the solution of the eigenvalue analysis. The frequency of this mode and consequently also the pressure distribution is highly influenced by the amount of heat release rate fluctuations considered in the model as shown in [START_REF] Laera | Impact of heat release distribution on the spinning modes of an annular combustor with multiple matrix burners[END_REF] (see the trajectory maps in Fig. 9 ). This results in an iterative procedure where the solution at the k th iteration is obtained using the ψ distribution computed at the k th -1 iteration. This procedure is iterated until the maximum error defined by

= | ψ ( ) k -ψ ( ) k -1 | is
lower than a threshold value, typically 10 -3 . It should be noticed that at each k th iteration, the nonlinear eigenvalue problem defined by Eq. ( 4) is solved in order to compute the new mode structure ψ [START_REF] Laera | Impact of heat release distribution on the spinning modes of an annular combustor with multiple matrix burners[END_REF] .

A validation of the proposed numerical procedure is discussed in Appendix A. The case investigated considers an annular cavity with a uniform distribution of heat release and a simplified model for the nonlinear flame response to pressure disturbances. Analytical solutions were derived for both spinning and standing limit cycles together with conditions for their stability [START_REF] Noiray | Investigation of azimuthal staging concepts in annular gas turbines[END_REF] . Confrontations between numerical simulations carried out in this work and analytical expressions perfectly match for both the rotating and standing modes validating the numerical procedure.

Experimental setup and numerical representation

We now briefly describe the MICCA combustor shown in Fig. 1 a. This system comprises an annular plenum connected by 16 injectors to an annular chamber formed by two cylindrical concentric quartz tubes of 200 mm length. Each injector consists of a cylinder of d br = 33 mm diameter and l br = 14 mm length exhausting gases through a l inj = 6 mm thick perforated plate featuring 89 holes of d p = 2 mm diameter located on a 3 mm square mesh. The system fed by a propane/air mixture allows the stabilization of a set of laminar conical flames above the 16 injectors as shown in Fig. 2 for one injector. Two operating conditions are investigated. Figure 2 a shows flames stabilized above a single matrix injector for a stoichiometric mixture φ= 1, while flames obtained for a slightly richer mixture at φ= 1.11 shown in Fig. 2 b have a higher longitudinal extension. These images were taken for stable operation in a single matrix injector setup. When the MICCA combustor is operated at condition A, a thermo-acoustic instability coupled to a spinning mode with stable limit cycle is observed [START_REF] Bourgouin | Characterization and modeling of a spinning thermoacoustic instability in an annular combustor equipped with multiple matrix injectors[END_REF] . When the MICCA is operated at condition B, a stable limit cycle coupled to a standing mode is found [START_REF] Durox | Nonlinear interactions in combustion instabilities coupled by azimuthal acoustic modes[END_REF] . Slanted [START_REF] Bourgouin | A new pattern of instability observed in an annular combustor: the slanted mode[END_REF] self-sustained combustion oscillations coupled to azimuthal modes were also identified in this setup when the flow operating conditions were modified.

The heat release rate is distributed in the simulations over a small volume located at the exit of each burner. It consists of a cylindrical volume of height l f and diameter d f . A previous study on the same combustor indicates that the dynamics of perturbations is influenced by the extension of the flame domain [START_REF] Laera | Impact of heat release distribution on the spinning modes of an annular combustor with multiple matrix burners[END_REF] . For conical flames, the flame volume dimensions are deduced by processing the images of the flame region under steady conditions shown in Fig. 2 as described in [START_REF] Laera | Impact of heat release distribution on the spinning modes of an annular combustor with multiple matrix burners[END_REF] . The numerical representation of the system is shown schematically in Fig. 3 a. The plenum consists of an annular cavity linked to the combustion chamber volume by sixteen injection units as in the real configuration. A model is used to represent the matrix injectors shown in Fig. 3 b. The body of each burner has the same dimensions as the real system. The perforated plate is replaced by a cylindrical volume having the same l inj = 6 mm thickness as the perforated plate and a base area with a diameter of d inj = 18.9 mm corresponding to the total flow passage area of the injector. The total height of the burner is l br + l in j = 20 mm. The combustion chamber is modeled as an annular duct open to the atmosphere with an augmented length of 41 mm to account for an end correction resulting in a total length of l cc = 241 mm. The value of this correction was determined experimentally by submitting the MICCA chamber to harmonic acoustic excitations near the unstable resonant mode and by scanning a microphone along a longitudinal axis above the combustion chamber annulus [START_REF] Laera | Impact of heat release distribution on the spinning modes of an annular combustor with multiple matrix burners[END_REF] . All other boundaries are treated as rigid adiabatic walls. The combustor operates at atmospheric conditions. The temperature of the plenum is set equal to 300 K.

Following previous studies [START_REF] Palies | Nonlinear combustion instability analysis based on the flame describing function applied to turbulent premixed swirling flames[END_REF][START_REF] Kinsler | Fundamentals of acoustics[END_REF] , the damping rate is deduced from a resonance response of the system by imposing an external perturbation with a loudspeaker and measuring the resonance sharpness. These measurements were carried out under cold flow conditions avoiding any corrections to account for absorption or generation of acoustic energy by the flame [START_REF] Palies | Nonlinear combustion instability analysis based on the flame describing function applied to turbulent premixed swirling flames[END_REF] . This, however, introduces some uncertainty since the value of the damping rate under hot conditions may differ from that estimated at room temperature. Figure 4 shows two acoustic response curves measured in the plenum ( Fig. 4 a) and in the combustion chamber ( Fig. 4 b) providing the resonance frequency bandwidth f at half-power. The damping rate δ appearing in Eq. ( 2) is deduced from the frequency bandwidth 2 δ = f ( s -1 ).

The azimuthal thermo-acoustic instabilities coupled by spinning and standing modes are discussed in what follows. Numerical simulations are compared with experiments and the stability of the numerical predictions is evaluated.

Analysis of operating point A

For operating condition A corresponding to a stoichiometric propane/air mixture φ= 1 with a bulk flow velocity measured in the cylindrical body of each injector equal to u b = 0.49 m s -1 , the system sustains a well-established spinning limit cycle coupled to the first azimuthal mode at a frequency of 487 Hz [START_REF] Bourgouin | Characterization and modeling of a spinning thermoacoustic instability in an annular combustor equipped with multiple matrix injectors[END_REF] .

In the numerical model, effects of steady combustion are taken into account through a temperature distribution in the gas stream. The temperature was measured, along a longitudinal axis passing in the center of one burner location, with a movable thermocouple and varies from 1470 K near the flame zone to 1130 K at the end of the combustion chamber. The flame domain consists of a cylindrical volume of height l f = 4 mm and diameter d f = 36 mm ( Fig. 3 b). These dimensions are deduced by processing the image shown in Fig. 2 a as in [START_REF] Laera | Impact of heat release distribution on the spinning modes of an annular combustor with multiple matrix burners[END_REF] . This gives a flame volume V f = 4.18 cm 3 which, for a global thermal power per burner of ˙ Q = 1.44 kW, yields a heat release rate per unit volume equal to ˙ q = 3 . 2 × 10 8 W m -3 . For each burner, the interaction between combustion and acoustics is expressed by making use of a global FDF determined experimentally in a single burner setup comprising the same injector and equipped with a driver unit, a hot wire and a photomultiplier to measure velocity and heat release rate fluctuations respectively [START_REF] Noiray | A unified framework for nonlinear combustion instability analysis based on the flame describing function[END_REF][START_REF] Boudy | Describing function analysis of limit cycles in a multiple flame combustor[END_REF] . The reference point for the velocity fluctuation measurements is located in the injection unit 20 mm below combustor backplane. Figure 5 shows the FDF used in this analysis. Measurements of the gain and phase lag were carried out for five velocity fluctuation levels ranging from | u / u | = 0.1 to | u / u | = 0.5 (white disks in Fig. 5 ). One difficulty is to gather FDF data at large forcing amplitudes. Due to limitations of the equipment used to modulate the flame, it was not possible to cover the full frequency and amplitude ranges.

A well-resolved FDF is used in the simulations by interpolation between the experimental points and extrapolation where experimental samples are missing. At very high amplitude levels that are reached in the MICCA experiment, the flames are disrupted and it is reasonable to represent this behavior by a fast drop in their response to external perturbation. Data are ex- trapolated at these levels with the best fit curve of a smooth fourth order polynomial function using data gathered at lower levels. This is carried out for each forcing frequency. There is a certain amount of uncertainty introduced by this process, but it is based on measured data points while most of the theoretical investigations are based on simplified representations. The nonlinear n -τ models [START_REF] Dowling | Nonlinear self-excited oscillations of a ducted flame[END_REF][START_REF] Li | Time domain simulations of nonlinear thermoacoustic behaviour in a simple combustor using a wave-based approach[END_REF][START_REF] Laera | A weakly nonlinear approach based on a distributed flame describing function to study the combustion dynamics of a fullscale lean-premixed swirled burner[END_REF] , simplified third order polynomials of heat release as a function of pressure [START_REF] Noiray | Investigation of azimuthal staging concepts in annular gas turbines[END_REF] or time invariant nonlinear representations of heat release rate as a function of pressure [START_REF] Ghirardo | Weakly nonlinear analysis of thermoacoustic instabilities in annular combustors[END_REF] do not feature all the complexity of the nonlinear flame response considered in the present study. Heat release rate fluctuations are assumed to be driven in the MICCA annular combustor by the fluctuating mass flow rates due to axial velocity perturbations through the corresponding injector [START_REF] Wolf | Using les to study reacting flows and instabilities in annular combustion chambers[END_REF][START_REF] Krebs | Thermoacoustic stability chart for high-intensity gas turbine combustion systems[END_REF] . This assumption is reasonable as the injectors are well separated in the configuration explored and there is no visible mutual interaction [START_REF] Worth | Modal dynamics of self-excited azimuthal instabilities in an annular combustion chamber[END_REF] . The reference points for the velocity fluctuations in the numerical domain are taken on the axis of the burner at the same distance from the flame domain as in the experiments.

Dynamics of an initially spinning mode

In a first stage, following the experimental observations the stability analysis is carried for the first azimuthal mode (1A) of the MICCA chamber starting the simulations with a spinning mode structure and the velocity distribution described by Eq. ( 6) . Figure 6 compares the system trajectories plotted in a frequencygrowth rate plane for two different values of the damping rate δ= 0 s -1 and δ= 12.5 s -1 in Eq. ( 2) . Each curve is colored with respect to the velocity fluctuation level that varies from | u / u | = 0.10 to | u / u | = 0.63 by steps of 0.01. Symbols are only plotted every 10 increments and at the trajectory endpoint to ease reading. The dynamical trajectories of the system are controlled by three free parameters, the frequency f , the growth rate α and the relative velocity amplitude level | u / u | . The 1A mode is found to be linearly unstable. For small velocity perturbations, the system features the highest growth rate of about 280 s -1 and a frequency around 400 Hz. A reduction of the growth rate and a substantial increase of the instability frequency is observed when | u / u | is augmented. This is due to the reduction of the FDF gain when the velocity fluctuation level increases, as can be observed in Fig. 5 a. If no damping is considered, the limit cycle condition is reached at a velocity fluctuation level that nullifies the gain of the FDF. When damping is considered the limit cycle is reached for an amplitude level | u / u | = 0.61, which corresponds to a FDF gain G = 0.08 at a frequency of f = 473 Hz. This frequency is close to that observed at limit cycle in the experiments f = 487 Hz. It is worth noting that at the velocity fluctuation level that nullifies the FDF gain, the sys- tem shows a negative growth rate of α= -12.5 s -1 which equals the damping rate that was determined experimentally. This checks that the dissipation rate is well represented in the numerical procedure.

The structure of the unstable mode is now investigated at the limit cycle. Since the sixteen burners have the same FDF and are assumed to operate at the same amplitude level, the circumferential symmetry of the system defined by the annular geometry is conserved. The nonlinear stability analysis leads to degenerate solutions featuring two azimuthal modes sharing the same frequency and spatial structures shifted by π /2 as in the linear case [START_REF] Bauerheim | Symmetry breaking of azimuthal thermo-acoustic modes in annular cavities: a theoretical study[END_REF] . It is thus possible to add these two solutions at each amplitude level and obtain a spinning mode.

The result is shown in Fig. 7 b (bottom) in the form of a pressure phase evolution plotted along the azimuthal direction at a mid-height position in the plenum and in the combustion chamber backplane. The phase evolutions feature a shift of 0.14 rad between the plenum and the combustion chamber. This angular shift is also observed in experiments and confirmed theoretically [START_REF] Bourgouin | Characterization and modeling of a spinning thermoacoustic instability in an annular combustor equipped with multiple matrix injectors[END_REF][START_REF] Laera | Impact of heat release distribution on the spinning modes of an annular combustor with multiple matrix burners[END_REF] . The pressure magnitude | ˆ p | shown in Fig. 7 a is obtained by plotting pressure contour lines computed by the Helmholtz solver on a cylindrical surface passing through the middle of the combustion chamber, the plenum, the burners and the microphone waveguides. The pressure iso-lines are deformed in the vicinity of the burners due to the near field acoustic interactions with the injectors, heat release zone and waveguides. However, a spectral analysis of the pressure distribution p = p (θ ) , not shown here, features small harmonic levels. At the burners locations, the harmonic content remains within 6% of the signal amplitude and falls to 1% one centimeter away from the chamber backplane as highlighted by the pressure iso-lines plotted in Fig. 7 a. This is confirmed observing the pressure distribu-

tion p (θ ) = | ˆ p | cos arg( ˆ p )
along the azimuthal direction plotted in Fig. 7 b (top). Deformations appear in the distribution taken at the backplane of the combustion chamber, whereas the influence vanishes in the plenum.

Stability of the limit cycle

The stability of the spinning equilibrium point is now investigated following ideas developed in [START_REF] Ghirardo | Weakly nonlinear analysis of thermoacoustic instabilities in annular combustors[END_REF] : limit cycles coupled to spinning modes are stable if the derivative of the FDF with respect to the amplitude of the spinning mode oscillation ( | u / u | sp ) is negative. For a constant damping rate, the derivative of the FDF around the equilibrium point can be approximated by the derivative of the growth rate, indicating that the stability criterion can be reformu-

∂α ∂| ˆ p | | u / u | sp < 0 (8)
The growth rate α plotted in Fig. 21 c as a function of the perturbation level u / | ū | features a negative derivative around the spinning oscillation equilibrium point reached for u / | ū | = 0.61. Eq. ( 8) is thus satisfied and one can conclude that the predicted spinning mode is stable as observed in the experiments.

Comparisons with experiments

Spinning pressure oscillations at the limit cycle corresponding One may however note that this harmonic content remains weak. The pressure peak at 974 Hz reaches 10 Pa. In the numerical calculation, signals are reconstructed by considering only the first harmonic found at 473 Hz in the simulations. Nevertheless, a good match is found in terms of amplitude and phase for the four sensors indicating that the fundamental dominates and proving that the numerical procedure is able to predict the difference in amplitude levels between the two cavities of the system. Again, a small phase mismatch can be observed when the comparison is carried out over a longer duration.

Time resolved heat release rate signals are deduced in the experiments from two photomultipliers (H1 and H2) equipped with an OH * filter and placed at locations shown in Fig. 1 c. The heat release rate fluctuation ˙ Q can also be deduced from the simulation:

˙ Q ˙ Q = ˆ u u F (ω r , | u / u | ) exp (iω r t ) , (9)
where ˆ u / ū is the calculated velocity oscillation level at limit cycle, F is the FDF that needs to be evaluated at the same forcing level, ˙ Q is the mean heat release rate, ω r is the angular frequency at limit cycle and t denotes the time. Figure 9 shows the recorded heat release rate signals (dashed lines) and the numerical reconstructions (solid lines). The two photomultipliers records feature nearly the same amplitudes and a phase shift of 1.63 rad being close to the theoretical value of π /2. However, the amplitude is

Table 1

Sensitivity of the limit cycle oscillation frequency f , the pressure amplitude recorded by microphone MC1 in the chamber, the FDF gain and phase lag and the relative level of heat release rate measured by photomultiplier H1 with respect to the velocity fluctuation level | u / ū | reached at limit cycle. continuous lines with rectangular marks in Fig. 9 . They feature the same phase shift as that recorded by the photomultipliers. In terms of amplitude, the two reconstructions share the same amplitude, but the value is about four times lower than the one recorded in the experiments. The reason for this sizable difference is now investigated with the help of Eq. ( 9) . Given the relatively low damping of the system ( δ = 12 . 5 s -1 ), the limit cycle is reached for a large velocity oscillation level | u / u | = 0.61 higher than 0.5, in a range where the FDF gain is slightly extrapolated and features a steep slope as can be seen in Fig. 5 1 . This is not the case for the corresponding heat release rate oscillation ˙ Q / ˙ Q due to the rapid drop of the FDF at large perturbation amplitudes. Data in Table 1 indicate that a reduction of 5% of the level of velocity fluctuation at the limit cycle results in a variation of 170% of the FDF gain G and, correspondingly, in the resulting heat release rate oscillation amplitude. The phase of the FDF and the frequency of the resonant mode remain as a first approximation unaffected by these changes. The amplitude differences between measurements and numerical predictions reduce when the heat release rate signals are reconstructed for an oscillation level | u / u | = 0.58 as shown by the continuous lines with circular marks in Fig. 9 . A perfect match in amplitude between the experimental and numerical signals is still not achieved, but differences are notably reduced. These tests confirm the strong sensitivity of the predicted level of heat release rate reached at limit cycle due to small uncertainties on the velocity oscillation level in the injector. This feature reflects that small uncertainties on the data gathered for FDF at high forcing levels may lead to large deviations of the predicted heat release rate fluctuations due to the rapid drop of the FDF gain with the forcing amplitude when the flame is disrupted.

| u / u | f (Hz) MC1 (Pa) Gain ϕ ( × π) | ˙ Q / ˙ Q | 0.

Analysis of operating point B

For an equivalence ratio φ= 1.11 and a bulk flow velocity equal to u b = 0.66 m s -1 , the system features a well-established selfsustained combustion oscillation associated to standing mode at a frequency of 498 Hz that was fully characterized in [START_REF] Durox | Nonlinear interactions in combustion instabilities coupled by azimuthal acoustic modes[END_REF] . Two standing mode patterns have been observed in the system for the same operating conditions depending on the run recorded. Figure 10 shows long time exposure photographs of these two modes with the position of the nodal line indicated by the red dashed line. The mode structure shown in Fig. 10 a features a nodal line between burners I -II and between burners IX -X . In this work, this mode structure will be designated as "V type". Figure 10 b shows the mode structure with the nodal line between burners V -VI and burners XIV -XV that will be referred as "H type".

The FDF corresponding to this new operating condition is shown in Fig. 11 . The same difficulties to get data for the FDF at high amplitudes, typically | u / ū | > 0 . 5 , persist for this new operating condition. As for operating condition A, the heat release rate is uniformly distributed over each burner by post-processing the flame image shown in Fig. 2 b [START_REF] Laera | Impact of heat release distribution on the spinning modes of an annular combustor with multiple matrix burners[END_REF] recorded for operating regime B in the single burner setup in a thermo-acoustically stable state. This procedure results in a flame volume V f = 6.11 cm 3 that is distributed over a cylinder of height l f = 6 mm, 2 mm longer than that used for operating regime A, and a diameter of d f = 36 mm which, for a global thermal power per burner of ˙ Q = 2.08 kW, yields a heat release rate per unit volume equal to ˙ q = 3.31 W m -3 . The mean temperature in the combustion chamber varies from 1521 K near the flame zone to 1200 K at the end of the combustion chamber. The temperature of the plenum is set equal to 300 K.

Dynamics of an initially standing mode

The numerical procedure described in Section 2 is used for the analysis of the 1A standing mode dynamics when the system operates at regime B. Without unsteady heat release, the circumferential symmetry of the system defined by the annular geometry is conserved and the eigenvalue analysis yields degenerate solutions. Figure 12 2 , computed under passive flame conditions ( ˙ q = 0) and plotted over a plane located at the burner inlet section. The two modes share the same frequency 472 Hz, but their structures are shifted by π /2. The corresponding distribution used to initialize the simulation is shown in Fig. 12 b where symbols indicate the burners' angular locations. When a distribution of FDFs is introduced ( ˙ q = 0), the rotational symmetry of the system is broken [START_REF] Noiray | Investigation of azimuthal staging concepts in annular gas turbines[END_REF] and the eigenvalue analysis results no more in degenerate azimuthal modes but yield two distinct waves characterized by different modal structures. The frequency and the growth rate of these modes depend on the level of asymmetry considered in the system.

At the onset of instability, Fig. 11 (b) indicates that the velocity fluctuation level weakly influences the flame response. The eigenvalue analysis of this weakly asymmetric system yields two waves with different yet close frequencies and growth rates (see Appendix B). Increasing the velocity fluctuation level, the differences increase between the gain and phase values taken by the FDFs from the different burners. This leads to stronger asymmetric configurations in which the two modes resulting from the eigenvalue analysis are characterized by an important shift in growth rate and also frequency. The frequency shift depends on the amount of heat release rate fluctuation considered in the system [START_REF] Silva | Combining a Helmholtz solver with the flame describing function to assess combustion instability in a premixed swirled combustor[END_REF] . In the validation calculations discussed in the Appendix B with a simplified heat release response, the two waves manifest only a growth rate shift, while they share the same frequency. Here the two solutions feature both a frequency and growth rate shifts. As discussed in Section 2 , for each velocity level | u / u | j , the mode structure ψ k computed at the k th iteration is used in Eq. ( 7) to modulate the FDFs. However, in order to be consistent, only the mode structure closer to the one chosen at the first iteration (the ψ 0 1 or ψ 0 2 shown in Fig. 12 ) is followed until the convergence criterion described in Section 2 is satisfied.

The dynamical trajectories of the system are first tracked without considering any damping by setting δ= 0 s -1 in Eq. ( 2) to analyze the influence of the chosen initial distribution. These trajectories are plotted in Fig. 13 the different levels, symbols are only plotted every 10 increments and at the trajectory endpoint. Rectangular marks in Fig. 13 indicate results obtained by considering the initial mode structure ψ (0) 1 shown in Fig. 12 a and designated as H type. Circular marks in Fig. 13 indicate results obtained by considering the initial mode structure ψ (0) 2 in Fig. 12 a and designated as V type. Regardless of the eigenmode used to initialize the simulations the frequencies and growth rates are the same. It is worth noting that, in contrast to spinning mode calculations, a limit cycle condition cannot be reached without accounting for a finite damping level regardless of the velocity oscillation level considered in the system. At each point of the trajectory, flames close to the nodal line always experience a small velocity fluctuation | u / u | . This small oscillation leads to a finite heat release rate fluctuation and this causes the mode to become unstable.

With the introduction of a damping term, a limit cycle is obtained as shown in Fig. 14 where the dynamical trajectories without damping, δ= 0 s -1 (circular marks), and with a damping rate of δ= 12.5 s -1 (rectangular marks) are compared. Each curve is colored in this figure by the velocity fluctuation level that varies from | u / u | = 0.1 to | u / u | = 0.9 with steps of 0.01. The limit cycle condition α= 0 is reached in the simulation for a fluctuation level | u / u | = 0.86 at a frequency f = 478 Hz. This value is close to that recorded at the limit cycle in the experiments f = 498 Hz. As discussed in [START_REF] Noiray | Investigation of azimuthal staging concepts in annular gas turbines[END_REF] , the amplitude of a standing limit cycle is greater than the amplitude of a spinning limit cycle that would settle for the same operating conditions and the same flame response model. Even though the operating conditions differ, it is found here that the oscillation level | u / u | = 0 . 86 of the limit cycle of operat- ing condition B associated to standing mode is larger than the one | u / u | = 0 . 61 found for the operating condition A featuring a limit cycle with a spinning pattern.

In agreement with experiments shown in Fig. 10 , two different modal structures with a nodal line shifted by π /2 are predicted at the limit cycle depending on the chosen initial condition. Figure 15 a shows the H type limit cycle pressure distribution calculated by the Helmholtz solver. This solution is obtained from the modal pressure distribution ψ (0) 1 shown in Fig. 12 as initial condition. Figure 15 b shows the V type limit cycle pressure distribution obtained assuming the initial distribution ψ (0) 2 . The distribution of the velocity fluctuation level | u / u | reached at limit cycle for both modal structures is shown in Fig. 16 . It differs from a pure sinusoid due to the local deformations of the pressure field near the injector outlets as already discussed for operating condition A. These deformations are taken into account in the ψ function used to fix the velocity oscillation level in the different burners. Symbols in Fig. 16 indicate the levels reached at the sixteen burner positions. One may note that for the largest velocity fluctuations, flow reversal conditions can be reached during part of the oscillation cycle in the injectors located close to the pressure anti-nodal lines.

It is next interesting to compare the location of the pressure nodal line. There is no precise experimental determination of its angular position. Examining the long time records shown in Fig. 10 , one finds that in both configurations the nodal line is always located between two burners in a region indicated in grey in Fig. 17 . The same figure also shows the predicted nodal line plotted as a red dashed line for the two structures. For the standing limit cycle featuring a V type structure, the nodal line is predicted between two burners in Fig. 17

Stability of the limit cycle

The stability of the standing limit cycle is now investigated. According to Ghirardo et al. [START_REF] Ghirardo | Weakly nonlinear analysis of thermoacoustic instabilities in annular combustors[END_REF] , a standing limit cycle is stable if three necessary and sufficient conditions are satisfied. In essence, the first condition requires that at the limit cycle, the growth rate decreases when the velocity level is increased. The growth rate trajectory obtained by varying the maximum amplitude of the velocity oscillation level | u / u | plotted in Fig. 23 c (continuous curve with rectangular marks) indicates that the derivative of the growth rate around the standing equilibrium point is negative. This test confirms that the first stability condition is satisfied. The second criterion defines a condition on the orientation of the nodal line. Ghirardo et al. [START_REF] Ghirardo | Weakly nonlinear analysis of thermoacoustic instabilities in annular combustors[END_REF] have demonstrated that this condition is always fulfilled in configurations with a large number of burners. Although the MICCA combustor features a finite number of burners, this second criterion may be regarded as satisfied considering that heat release rate fluctuations take place in a flame domain, which nearly covers the entire surface area of the combustion chamber. The third condition discusses the stability of the standing wave pattern. Stable standing modes need to comply with the following inequality:

N 2 n = 2 π 0 F (ω r , | u / u | i ψ (θ )) /Z(θ ) sin (2 nθ ) dθ > 0 . ( 10
)
In this expression, n is the azimuthal mode order, which is here equal to n = 1 and the mode ψ is chosen such that there is a pressure anti-node at θ = π / 4 . The H type mode has the closest structure fulfilling this condition and is used for the calculation. The FDF values F in Eq. ( 10) are calculated at the limit cycle oscillation frequency f = 480 Hz with the corresponding velocity distribution.

The quantity Z ( θ ) designates the impedance, i.e., the ratio of pressure in the flame zone to the velocity at the reference point.

Results for

F (| u / u | i ψ (θ )) /Z(θ )
are plotted in red in Fig. 18 .

It consists of a piecewise function taking constant values over the angular extensions of the flame zones and zero values in the angular extensions between two flames. In the same graph the function sin (2 n θ ) is represented in black. The consequence is that the component in red in Fig. 18 changes sign with amplitude resulting in a positive overall integral for N 2 n . The predicted standing limit cycle is thus found to be stable as observed in experiments. This is due to the fact that the FDF of the matrix burners investigated in the present study features a phase lag which is sensitive to | u / u | as shown in Fig. 11 b. The situation slightly differs from that investigated in Appendix A with a simplified flame response characterized by a monotonically decreasing gain but a constant phase lag independent of the input amplitude level. In this case, limit cycles coupled to standing waves are found to be unstable.

Comparisons with experiments

Calculations for the standing limit cycle are now compared to measurements. The frequency and amplitude of velocity fluctuations predicted at limit cycle at one reference point of the system are used to reconstruct time resolved pressure fluctuations in the plenum and in the chamber and heat release rate fluctuations. The analysis is only carried out for results corresponding to the V structure type. Figure 19 (a The amplitude and the phase shift between microphones located near the anti-nodal line, i.e., MP1 and MP5, is well captured by the simulations. Calculations also retrieve the amplitudes of the microphones located near the nodal line, i.e., MP3 and MP7, however, with a small phase mismatch. This is due to the shift of θ ∼ 3 °between the angular position of the sensors and of the nodal line shown in Fig. 17 a. In the numerical reconstruction, microphones MP3 and MP7 are located on the two sides of the nodal line and their signals are consequently in phase opposition. In the experiments, these records are nearly in phase as shown by the green and blue dashed lines in Fig. 19 (a) (top). Increasing the interrogation period one finds a small phase mismatch between experiments and simulations due to the 20 Hz difference between the measured and predicted limit cycle oscillation frequencies.

Figure 19 (a-bottom) shows the same type of comparison for pressure signals in the combustion chamber. Again, numerical reconstructions and measurements are, respectively, identified by continuous and dashed lines. In the experiments, the signals feature a some harmonic content and significantly differ from one to another depending on the azimuthal position of the sensors. The analysis of the spectrum content of the signals shown in Fig. 19 (bbottom) reveals the second harmonic peak at 992 Hz. A previous study on the same combustor featuring standing mode oscillations [START_REF] Durox | Nonlinear interactions in combustion instabilities coupled by azimuthal acoustic modes[END_REF] has shown that this frequency is associated with the first longitudinal mode of the combustion chamber (1L). For the pressure signals located close to the nodal line, the amplitude of the second harmonic peak is comparable to the amplitude of the 1A mode putting in evidence a competition between the 0L-1A plenum oscillation at 495 Hz and a longitudinal mode associated to the 1L-0A mode of the chamber. This yields the distorted signals recorded by microphones MC3 and MC7 shown in Fig. 19 (a-bottom). The other microphones located closer to the pressure anti-nodes feature a second harmonic at 992 Hz with an amplitude of about one order of magnitude lower than that of the first harmonic [START_REF] Durox | Nonlinear interactions in combustion instabilities coupled by azimuthal acoustic modes[END_REF] . This yields the roughly sinusoidal signals with a peak amplitude at 495 Hz of 60 Pa recorded by microphones MC1 and MC5 shown in Fig. 19 (a-bottom). In the FDF framework, the numerical signals can only be reconstructed for the first harmonic oscillation. Overall, a good match is found in terms of both amplitude and phase for the four sensors indicating that the numerical procedure is able to pre- dict the difference in amplitude levels between the two cavities of the system. Again, a small phase mismatch is observed when the comparison is carried out over a longer duration period.

The numerical heat release rate signals (continuous lines) reconstructed with Eq. ( 9) are compared in Fig. 20 with the two OH * light intensity signals (dashed lines) recorded by the photomultipliers (see Fig. 1 (c)). The two photomultipliers signals feature different amplitudes. The injector close to the nodal line features the largest heat release rate oscillations reaching ˙ Q / ˙ Q = 0 . 5 measured by H1. A drop of about 50% of the peak amplitude is observed for the injector close to the anti-nodal line (H2 signal in Fig. 20 ). It is also worth noting that the signal is in this case not a pure sinusoidal wave indicating the presence of harmonics. The two signals measured by the photomultipliers are nearly in phase due to the position of the sensors with respect to the nodal line. The small phase shift is due to the presence of harmonics. As for the pressure signals, the simulated heat release rate signals only consider the first harmonic. The reconstructed numerical signals have the same phase shift as those recorded by the photomultipliers. This indicates that the phase of the FDF is well captured. The numerical procedure also retrieves the amplitude drop between the two signals, however, the predicted amplitude peak and the one recorded in the experiments differ. For the H1 signal, differences between simulations and experiments are mainly due to the mismatch on the position of the nodal line shown in Fig. 17 . As already discussed for operating condition A, the heat release rate oscillation level strongly depends on the predicted value of | u / u | and on the gain of the FDF in Eq. ( 9) . For high velocity fluctuation levels, the FDF gain rapidly drops as shown in Fig. 11 a. As a consequence, the numerical prediction of the H2 signal is quite sensitive to small uncertainties on the velocity oscillation level at the limit cycle. These uncertainties lead to important variations of the FDF gain G and consequently in the predicted heat release rate oscillation at limit cycle. The origin of the differences observed at limit cycles between measurements and simulations is thus the same as for the spinning mode calculations studied in the previous section.

Mode type selection

Following the suggestion of a reviewer, an investigation of the possible scenarios leading to the spinning and standing limit cycles analyzed in the previous sections is now proposed.

Analysis for the operating point A

The stability analysis is repeated for operating condition A, but instead of assuming an initially spinning mode as done in Section 4.1 , simulations are initiated with a standing mode structure and the velocity distribution described by Eq. ( 7) . The corresponding trajectory (black continuous line with rectangular marks) plotted in a frequency-growth rate-velocity fluctuation level space is compared in Fig. 21 a with the spinning mode trajectory calculated in Section 4.1 (red continuous line with circular marks). At the starting point, i.e. for | u / u | = 0.1, the two trajectories perfectly match. This is due to the fact that in the linear regime, i.e. for values of | u / u | ≤0.1, the FDF is simply a transfer function and its gain and phase do not depend on the input amplitude level. This also means that in simulations started with a standing mode structure, each burner of the chamber operates with the same velocity fluctuation. Projecting these trajectories on a frequency-growth rate ( f -α) plane, as shown in At the equilibrium point one finds that the two solutions share the same frequency, but differ in the level of velocity fluctuations: the standing mode limit cycle is obtained at higher values of | u / u | with respect to the spinning mode. One may then first discuss the stability of the two equilibrium points. The spinning limit cycle has been already proved to be stable in Section 4.2 . To prove the stability of the standing mode solution, the three criteria summarized in Section 5.2 are applied to this new standing equilibrium point. This analysis with the same experimental FDF shown in Fig. 5 has been recently carried out by Ghirardo et al. [START_REF] Ghirardo | Weakly nonlinear analysis of thermoacoustic instabilities in annular combustors[END_REF] . According to this study, a standing limit cycle assumed to occur for a value of | u / u | st = 0.5 and for a frequency of f = 487 Hz has been shown to be unstable and this was accepted as the main reason why the spinning mode prevails in this case. One may note, however, that the amplitude level adopted in [START_REF] Ghirardo | Weakly nonlinear analysis of thermoacoustic instabilities in annular combustors[END_REF] differs from the value obtained from the present numerical results where | u / u | st 0 . 86 and the frequency is slightly at variance with that deduced from the present calculations. The stability analysis of the limit cycle at | u / u | st 0 . 86 and f = 473 Hz is thus repeated here. The standing mode remains unstable as observed from Fig. 22 where inequality of Eq. ( 10) is graphically verified as described in Section 5.2 .

In addition to the stability analysis of the equilibrium solutions, a second possible mode type selection scenario may be considered by analyzing the two trajectories shown in Fig. 21 c. One may observe that the growth rate of the standing mode trajectory is always greater than the one of the spinning mode, with the exception of the range with | u / u | → 0.1 where the spinning mode trajectory has a greater growth rate. This region is highlighted in the zoomed box in Fig. 21 c. Noting that the two dynamical trajectories do not cross ( Fig. 21 a), one may conclude that the nature of the limit cycle is determined in the region where the FDF begins to depend on the input amplitude. According to this scenario, the spinning mode of oscillation is selected in this case because its growth rate slightly exceeds that of the standing mode. Once on that trajectory the velocity fluctuation level remains uniformly the same for all injectors and this leads to a spinning mode limit cycle. Increasing the value of | u / u | on the spinning mode path all the flames respond to the same velocity fluctuation level. A jump to the standing mode solution might be possible in principle but this would require a strong perturbation that would change the velocity distribution and provide the nonuniformity in amplitude pertaining to the standing mode oscillation. The difference in growth rates of the standing and spinning solutions is relatively small and one cannot be totally sure that this explains why the spinning mode is experimentally observed in this case. One can say at least that the range where the bifurcation takes place corresponds to the region where the FDF map shown in Fig. 5 is interpolated between experimental data. The errors in that region are lower than those made when the amplitude of oscillation becomes large near the limit cycle of the standing mode.

One may conclude from this analysis that the selection between spinning and standing modes may be explained by examining the stability of the final limit cycles in agreement with criteria derived by Ghirardo et al. [START_REF] Ghirardo | Weakly nonlinear analysis of thermoacoustic instabilities in annular combustors[END_REF] . However, the mode type selection is perhaps not only a matter of limit cycle stability. Simulations also reveal that the spinning mode features a growth rate exceeding that of the standing mode in a narrow region where the FDF starts to change with the input amplitude level. While the difference is minute, an alternative scenario has been developed to explain the spinning mode selection. A jump between the spinning and standing limit cycles seems to be possible but would require a sufficiently strong perturbation modifying the azimuthal distribution in velocity fluctuation. This mechanism would then be somewhat similar to the mode switching process documented in Noiray et al. [START_REF] Noiray | A unified framework for nonlinear combustion instability analysis based on the flame describing function[END_REF] .

Analysis for operating point B

As for the previous configuration, the stability analysis is now repeated for operating point B with an initial spinning mode structure and the velocity distribution in Eq. ( 6) . The predicted trajectory (red continuous line with circular marks) is compared with the trajectory calculated in Section 5.1 in the three-dimensional 23 c, one finds that the growth rate deduced from the standing mode distribution is always greater than the growth rate of the spinning mode distribution. According to this criterion, the standing mode prevails in this situation leading to a limit cycle with a standing pattern in agreement with experiments.

f -α -| u / u | space in
The stability of the two possible equilibrium points is now assessed. The standing equilibrium point was already proved to be stable in Section 5.2 . However according to Eq. ( 8) , the spinning equilibrium point is also found to be stable. This state is however not observed in the experiment. As a conclusion, the observed stable limit cycle state is in agreement with the stability criteria developed in Ghirardo et al. [START_REF] Ghirardo | Weakly nonlinear analysis of thermoacoustic instabilities in annular combustors[END_REF] , but these criteria are insufficient for this operating condition to fully determine which state is observed in the experiment. This case indicates without ambiguity that the dynamical trajectories need to be considered to fully determine the state of oscillation observed in the experiments. The main findings for the stability and trajectory analyses of the two operating points are summarized in Table 2 .

It is concluded that a criterion combining the stability analysis at the limit cycle with the trajectory analysis (line 7 in Table 2 ) might best define the mode type at the limit cycle.

Conclusion

A procedure combining a Helmholtz solver and the FDF framework is proposed to analyze the dynamics of an annular combustor. The FDF is determined experimentally, together with the damping rate. This method is first validated by making use of an idealized annular combustor model characterized by a simplified nonlinear flame response. Subsequently, the dynamics of spinning and standing combustion instabilities observed in a laboratory scaled annular combustor are then investigated numerically. Spinning and standing initial solutions are considered first without effects of unsteady heat release. For the spinning mode dynamics, effects of unsteady heat release is treated by assuming that the velocity oscillation level | u / u | is the same in each burner and that the flame response described by the FDF operates at this level in the 16 injectors.

For the standing mode dynamics, all injectors operate at different velocity oscillation levels | u / u | (θ ) depending on their relative position θ with respect to the nodal line of the analyzed mode. The distribution of | u / u | over the injectors is deduced from the pressure distribution computed by the code using an iterative procedure. Calculations of the dynamical trajectories of spinning and standing modes are carried out for increasing velocity oscillation levels until a limit cycle is reached when the calculated growth rate equals the damping rate. This numerical procedure is used to tentatively determine the modal structure at the limit cycle. Results indicate that the system trajectories and the stability of the limit cycle need both to be considered to fully determine the final state observed in the experiments for the two operating conditions explored. The problem of mode type selection however deserves further investigation, since the present results correspond to special cases with no proof of generality.

Results at limit cycle for pressure, velocity and heat release rate signals are for the first time compared to detailed experimental data.

The main findings of this study are:

• The predicted instability frequencies match those measured in the experiments for the two different operating conditions investigated. • The proposed numerical strategy allows to capture the correct pressure signals recorded in the plenum and in the chamber with a good match in terms of phase shift and amplitude for both operating conditions featuring a stable spinning limit cycle and a stable standing limit cycle. • For the standing limit cycle, the pressure signals close to the nodal line are less well reproduced than those at the anti-nodal line, but the two possible positions of the nodal lines observed in the experiments between two burners are well captured. • The stability of the predicted spinning and standing limit cycles has been proven numerically by making use of recent theoretical elements and match the states observed in the experiments. • Heat release rate signals at spinning and standing limit cycles are well reproduced in terms of phase shift and of heat release distribution between the different burners, but the predicted amplitudes differ from those measured. This is mainly attributed to the fact that the limit cycle is reached for oscillation states where the gain of the FDF rapidly drops and nonlinearly with the velocity fluctuation level | u / u | . This means that even small uncertainties on | u / u | lead to large variations of the FDF gain and, consequently, to large uncertainties on the amplitude of the predicted heat release rate fluctuations. • Two scenarios are investigated to deduce the preferred modal type that will settle at limit cycle. The first scenario based on a stability analysis of the predicted limit cycle oscillation distinguishes the final outcome in one case but fails to do so in another case. A new scenario relying on a calculations of frequency-growth rate trajectories of initially spinning and standing modes appears to identify the solution type. In one case the difference in growth rate is relatively small and may not be completely meaningful. The combination of these two scenarios might perhaps give the answer but this admittedly needs to be confirmed by further investigations.

This study shows that the numerical procedure developed herein adequately reproduces the dynamics of combustion instabilities coupled to azimuthal modes provided that the FDF can be determined accurately at sufficiently large perturbation levels. The lack of such data is a source of uncertainty.

Acknowledgment

The development of the annular combustor MICCA was funded by the Agence Nationale de la Recherche , Contract No. ANR-08-BLAN-0027-01 , and by Safran Tech. Davide Laera was supported by a research fellow grant provided by the Politecnico di Bari (Bari, Italy) for a six months period at EM2C. The authors wish to thank the reviewers for their critical and constructive examination of this article.

Appendix A. Validation of the numerical procedure

Theoretical results from Noiray et al. [START_REF] Noiray | Investigation of azimuthal staging concepts in annular gas turbines[END_REF] are used validate the proposed numerical procedure. The main elements of this analytical model are reproduced here. The geometry considered is an idealized annular combustion chamber of radius R where pressure fluctuations are assumed to be only function of time and azimuthal coordinate θ . After proper normalization of the flow variables, the normalized pressure disturbances p obey to:

∂ 2 p ∂t 2 + ζ ∂ p ∂t - ∂ 2 p ∂θ 2 = ∂Q ∂t , ( 11
)
where Q represents the normalized volumetric heat release rate disturbance and ζ accounts for damping. It is further assumed that the pressure of any azimuthal mode of order n can be expressed as:

p(θ , t ) = A (t ) sin (nt ) cos (nθ ) + B(t ) cos (nt ) sin (nθ ) , (12)
where A and B are slowly varying functions with time. For a nonlinear flame response linked to the pressure by:

Q = β p -κ p 3 , (13)
Noiray et al. [START_REF] Noiray | Investigation of azimuthal staging concepts in annular gas turbines[END_REF] found that the dynamics of this system was determined by a set of two-nonlinear coupled differential equations:

dA dt = (β -ζ ) A - κ 32 9 A 2 + 3 B 2 A , ( 14
)
dB dt = (β -ζ ) B - κ 32 9 B 2 + 3 A 2 B. ( 15
)
The spinning and standing limit cycles correspond to the fixed points of this system. Spinning modes are found for

A = B = 2 (β -ζ ) / 3 κ.
Standing modes correspond to solutions with one Figure 24 shows the computational domain used for the Helmholtz solver simulations. It consists of an annular duct with a length of 0.2 m. The radial extension of the annulus is set to 1 mm to avoid any low frequency radial components. Furthermore, a velocity node is assumed at each boundary to exclude the rising of longitudinal components. Under passive flame conditions ( Q = 0 ), the eigenmode analysis leads to degenerate solutions due to the circumferential symmetry of the system. In this case, the two first azimuthal modes ψ 1 and ψ 2 share the same frequency and are shifted by π /2. These two modes are orthogonal and form a basis used to describe the pressure field in the system. The heat release rate is assumed to be distributed in a flame sheet volume extending 4 mm in the longitudinal direction and covering the entire inlet section indicated in blue in Fig. 24 . The Helmholtz solver requires the transposition of the nonlinear flame model in the frequency domain. Taking the Fourier transform of Eq. ( 13) yields:

ˆ Q ω r , | ˆ p | = ˆ p ( ω ) 3 4 κ| ˆ p | 2 -β , ( 16
)
where ˆ Q and ˆ p are the dimensionless Fourier transforms of heat release rate and pressure disturbances taken at the inlet.

The weakly nonlinear stability analysis around the first azimuthal mode of the system (1A) is repeated for different pres- For the spinning mode calculation, the same forcing level is imposed at each injector. The circumferential symmetry of the system defined by the annular geometry is thus conserved. Consequently, regardless of the pressure fluctuation level | ˆ p | considered in the model, the weakly nonlinear stability analysis always yield degenerate solutions. This indicates that the two degenerate modes ψ 1 and ψ 2 have the same amplitude, a result in line with [START_REF] Noiray | Investigation of azimuthal staging concepts in annular gas turbines[END_REF] . also equal to that found analytically [START_REF] Noiray | Investigation of azimuthal staging concepts in annular gas turbines[END_REF] . If the pressure fluctuation level is increased further, damping of acoustic energy becomes greater than the energy gained and the oscillation amplitude ceases to grow. For | ˆ p | = 1 the FDF saturates at 3 / 4 κβ, the growth rate is negative (damping rate) and is equal to α = -ζ / 2 .

For the standing mode calculation, the value of the FDF gain depends on the angular position and is fixed by the distribution analytical predictions [START_REF] Noiray | Investigation of azimuthal staging concepts in annular gas turbines[END_REF] . The stability of the predicted limit cycles is now investigated. By inspection, one observes that the stability condition of Eq. ( 8) is satisfied for the growth rate behavior shown in Fig. 25 a, indicating that the simulated spinning limit cycle is found to be stable in agreement with analytical results [START_REF] Noiray | Investigation of azimuthal staging concepts in annular gas turbines[END_REF] . The stability of the standing mode is analyzed by determining the Jacobian matrix. Since B= 0 when A = 0 and vice versa, the only non-zero components of the Jacobian matrix are the coefficients: For each point of the red trajectory B= 0, indicating that this mode has no influence on the standing equilibrium point of this system [START_REF] Noiray | Investigation of azimuthal staging concepts in annular gas turbines[END_REF] . For a given oscillation level | ˆ p | max , the growth rate α ψ 2 value on the red curve corresponds to the Jacobian coefficient J 22 . At the equilibrium point J 22 is found to be equal to λ 2 = ( βζ ) / 3 in agreement with the analytical value [START_REF] Noiray | Investigation of azimuthal staging concepts in annular gas turbines[END_REF] . The eigenvalues λ 1 and λ 2 of the computed Jacobian matrix have opposite signs indicating that the standing mode is a saddle point and corresponds thus to an unstable limit cycle in agreement with theoretical predictions.

These simulations show that the proposed numerical procedure suitably retrieves the amplitude and the stability properties of both standing and spinning limit cycle in an idealized configuration where analytical results are available.

  plitudes of velocity fluctuations | u / u | have to be considered to represent the response of the flame above each injector. This leads to consider different FDF gains and phase lags values for each injector. In the proposed methodology, for each level of velocity fluctuations | u / u | j , the distribution of | u / u | over the injectors is determined from the pressure distribution computed by the code by setting:

Fig. 1 .

 1 Fig. 1. (a) Photograph of the MICCA combustor with a close up view of a matrix injector and a waveguide outlet. (b) Schematic representation of the experimental setup. (c) Top view of the MICCA chamber with microphone and photomultiplier measurement locations indicated.

Fig. 2 .

 2 Fig. 2. Images of the flame region recorded above one matrix injector under stable operation for the two conditions investigated.

Fig. 3 .

 3 Fig. 3. (a) Top view schematic representation of the model. (b) Longitudinal A-A cut showing geometrical details of the matrix injector model and of the flame domain. (c) Three dimensional model of the MICCA chamber with the details of the unstructured mesh comprising approx. 130,0 0 0 tetrahedral elements.

Fig. 4 .

 4 Fig. 4. Acoustic response of the MICCA combustor from 300 Hz to 500 Hz measured by microphones located in the (a) plenum and (b) combustion chamber. The frequency bandwidth f determined at half maximum provides the damping rate in both volumes.

Fig. 5 .

 5 Fig. 5. Interpolated Flame Describing Function (FDF) obtained for operating point A: φ= 1.00 and u b = 0.49 m s -1 . Experimental data points are displayed as white dots. (a) Gain. (b) Phase ϕ.

Fig. 6 .

 6 Fig. 6. Dynamical trajectories in the frequency ( f ) -growth rate ( α) plane colored by the velocity fluctuation level | u / u | for two damping rates δ = 0 (rectangular marks) and δ= 12.5 s -1 (circular marks). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

Fig. 7 .

 7 Fig. 7. (a) Pressure mode magnitude | ˆ p | with pressure contour lines plotted on a cylindrical surface equidistantly located from the lateral walls. (b) Top: pressure structure ˆ p = | ˆ p | cos arg( ˆ p ) distribution along the azimuthal direction in the plenum (circular marks) and in the combustion chamber (rectangular marks). Bottom: pressure phase evolution in the azimuthal direction in the plenum (circular marks) and in the combustion chamber (rectangular marks).

Fig. 8 .

 8 Fig. 8. (a) Four pressure signals recorded by microphones in the plenum (top) and in the combustion chamber (bottom) under spinning limit cycle conditions (dashed lines) compared with numerical reconstructions (solid lines). (b) Spectral content of the pressure signals in the plenum (top) and in the combustion chamber (bottom). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

Fig. 9 .

 9 Fig. 9. Time resolved heat release rate signals recorded by two photomultipliers in the combustion chamber (dashed lines) compared with numerical reconstructions (solid lines) for | u / u | = 0.61 (rectangular marks) and for | u / u | = 0.58 (circular marks). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

  to | u / u | = 0.61 are now compared to measurements. In the experiments, these fluctuations are recorded by microphones placed at four positions equidistantly separated on the external perimeter of the plenum and at four positions in the backplane of the combustion chamber as shown in Fig.1 c. Figure8(a) displays microphone measurements in the plenum (dashed lines) compared with the numerical reconstructions (continuous lines) at the same locations. The phase shift between two microphones signals corresponds to their relative position confirming the spinning nature of the mode. Well-established sinusoidal signals with a peak of about 260 Pa at a frequency of 487 Hz, as shown in Fig.8(b) (top), are found in the experiments. The amplitude and the phase shift between microphones is well captured in the simulation. Increasing the examination period, a small phase mismatch between experiments and simulations appears due to the 15 Hz difference between the numerical and experimental frequencies.

Figure 8 (

 8 bottom) shows experimental (dashed lines) and numerical signals (solid lines) in the combustion chamber. The microphones mounted on waveguides at a distance of 170 mm away from the backplane of the combustion chamber measure a delayed signal with a time lag τ m -b = 0.5 ms. Since this delay is not negligible compared to the oscillation period of the instability ( 2 ms), it is taken into account in the data processing. It is first worth noting that the pressure fluctuation level only reaches 60 Pa near the chamber backplane and experimental signals are not purely symmetric relative to the ambient mean pressure, indicating the presence of harmonics as revealed in the spectral content of the pressure signals shown in Fig. 8 (b) (bottom).

  not constant with time indicating the presence of harmonics. As for the pressure signals, the simulated heat release rate signals only consider the oscillation at the fundamental frequency. The reconstructed numerical signals with | u / u | = 0.61 correspond to the

Fig. 10 .

 10 Fig. 10. Long time exposure photographs of the flames at a limit cycle coupled to a stable standing azimuthal mode. The velocity nodal line is shown as a dashed red line. (a) V type mode structure with a nodal line observed between burners I -II and IX -X . (b) H type mode structure with with a nodal line between burners VIII -IX and XIV -XV [37] . (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

  a. A small change of the velocity fluctuation level | u / u | weakly alters the corresponding pressure oscillation level in the plenum and the chamber because these quantities scale linearly with | u / u | as shown in Table

Fig. 11 .

 11 Fig. 11. Interpolated Flame Describing Function (FDF) obtained for operating point B: φ= 1.11 and u b = 0.66 m s -1 . Experimental data points are displayed as white dots. (a) Gain. (b) Phase ϕ.

Fig. 12 .

 12 Fig. 12. (a) Pressure modulus | p | of the f = 472 Hz degenerate modes computed under passive flame conditions and plotted in a plane located at the burner inlet section. (b) Pressure distribution ψ used to initialize the computation of the standing mode. Symbols indicate the angular position of the sixteen burners.

  a shows the pressure modulus | p | of the degenerates modes, ψ (0) 1 and ψ (0)

  in a frequency-growth rate plane and colored with respect to the maximum amplitude level | u / u | j , which is varied from | u / u | = 0.1 to | u / u | = 0.9. The subscript j will be omitted in what follows. Unless indicated otherwise, | u / u | always refers to the maximum velocity fluctuation level considered in the simulation for the standing mode analysis. In order to emphasize

Fig. 13 .

 13 Fig. 13. Dynamical trajectories colored by the maximum velocity oscillation level | u / ū | in absence of damping δ= 0 s -1 . Square marks indicate results obtained for the H type standing mode structure. Circular marks indicate results obtained for the V type mode structure. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

Fig. 14 .

 14 Fig. 14. Dynamical trajectories in the frequency ( f ) -growth rate ( α) plane colored by the velocity fluctuation level | u / u | without damping δ= 0 (circular marks) and with damping δ= 12.5 s -1 (square marks) distributed uniformly in the numer- ical domain. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

Fig. 15 .

 15 Fig. 15. (a) Pressure modal distribution | ˆ p | / | ˆ p | max featuring a H type (a) or V type (b) mode structure calculated by the Helmholtz solver with pressure contour lines plotted on a cylindrical surface passing through the middle of the combustion chamber, the plenum, the burners and the microphone waveguides.

Fig. 16 .

 16 Fig. 16. Velocity distribution at the burners' positions for the standing mode at limit cycle with the V type (black dotted line) and H type (red dotted line) mode structure. Circular symbols indicate the burners' positions oriented as shown in Fig. 17 . (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

Fig. 17 .

 17 Fig. 17. Indication of the predicted angular position of the numerical pressure nodal line (red dashed line) together with the angular region (highlighted in grey) in which the experimental nodal line is recorded for both observed standing mode structures (a) V type (b) H type. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

  a. A shift of π /2 is observed for the H type structure in Fig. 17b.

Fig. 18 .

 18 Fig. 18. The piecewise function F (ω r , | u / u | i ψ (θ )) /Z(θ ) for the combustor under analysis (red line) together with a sin (2 n θ) (black line). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

  -top) displays pressure measurements in the plenum (dashed lines) compared with the numerical reconstructions (continuous lines) at the same locations in the numerical domain. Microphones MP3 and MP7, located at θ MP 3 = 78.8 °and θ MP 7 = 258.7 °, respectively, detect pressure oscillations close to the nodal line and they consequently feature low fluctuation levels in the experiments. Microphones MP1 and MP5 located near the anti-nodal line record a well-established sinusoidal signal with a peak amplitude of 350 Pa at a frequency of 498 Hz as shown in Fig. 19 (b-top). Given the standing nature of the recorded mode, the phase shift between two microphone signals corresponds to their relative position with respect to the nodal line. Records of microphones located at opposite sides of the nodal line are shifted by π .

Fig. 19 .

 19 Fig. 19. (a) Four pressure signals recorded by microphones in the plenum (a) and in the combustion chamber (b) at limit cycle (dashed lines) compared with numerical reconstructions (solid lines). (b) Spectral content of the pressure signals in the plenum (top) and in the combustion chamber (bottom). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

Fig. 20 .

 20 Fig. 20. Time resolved heat release rate signals for the standing mode recorded by H1 and H2 photomultipliers in the combustion chamber (dashed lines) compared with numerical reconstructions (solid lines).

Fig. 21 .

 21 Fig. 21. (a) Trajectories f -α -| u / u | of solutions of Eq. (2) for operating regime A initiated with a standing mode structure (black continuous line with rectangular marks) and a spinning mode structure (red continuous line with circular marks). Projection of the trajectories on f -α plane (b) and α -| u / u | plane (c). Circle markers indicate results for an initially spinning mode structure, square markers indicate results for a standing mode structure. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

Fig. 22 .

 22 Fig. 22. The piecewise function F (ω r , | u / u | i ψ (θ )) /Z(θ ) for the standing equilibrium point at | u / u | st 0 . 86 and f = 473 Hz (red line) together with a sin (2 n θ) (black line). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

  Fig. 23 a. As already observed during the analysis of the trajectories of the system dynamics operating in regime A, the two trajectories obtained in regime B perfectly match for low velocity fluctuation levels when | u / u | = 0.1. The paths then diverge leading to different limit cycles. In distinction with regime A studied in Fig. 21 b, the two trajectories plotted in the frequencygrowth rate plane in Fig. 23 b overlap only in the linear range. When | u / u | > 0 . 1 , the two trajectories diverge leading to limit cycles with different frequencies and amplitudes. Analyzing the trajectories in the α -| u / u | plane in Fig.

Fig. 23 .

 23 Fig. 23. Trajectories f -α -| u / u | of solutions of Eq. (2) for operating regime B initiated with a standing mode structure (black continuous line with rectangular marks) and a spinning mode structure (red continuous line with circular marks). Projection of the trajectories on f -α plane (b) and α -| u / u | plane (c). Circle markers indicate results for an initially spinning mode structure, square markers indicate results for a standing mode structure. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

  Prediction of trajectory analysis Spinning mode at f = 473 Hz with peak pressure of 260 Pa in the plenum and 57 Pa in the chamber Standing mode at f = 478 Hz with a peak pressure of 345 Pa in the plenum and 48 Pa in the chamber Criterion combining stability at limit cycle and trajectory analysis Spinning mode Standing mode Experimental limit cycle oscillation Spinning mode at f = 487 Hz with a peak pressure of 260 Pa in the plenum and 60 Pa in the chamber Standing mode at f = 498 Hz with a peak pressure of 350 Pa in the plenum and 60 Pa in the chamber

Fig. 24 .

 24 Fig. 24. Computational domain used for the validation of the numerical procedure. The flame domain is highlighted in blue. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

•

  sure amplitude levels starting from | ˆ p | = 0 and incrementing this value until a limit cycle condition is reached. Considering that the heat release rate is continuously distributed in the flame domain, and designating the angular coordinate by θ , one gets the following pressure distributions: Spinning mode: | ˆ p | ( θ ) =C, where C is a real positive constant, • Standing mode: | ˆ p | ( θ ) = | ˆ p | j ψ ( θ ) .

Figure 25 a

 25 shows results for the growth rate α as a function of the pressure fluctuation level | ˆ p | , which is varied from 0 to 1. For | ˆ p | → 0 the FDF corresponds to a linear flame transfer function (FTF). The value of the growth rate α = 1 / 2(β -ζ ) derived from the linear stability analysis from [16] is retrieved by the simulation when | ˆ p | → 0 . Increasing the pressure fluctuation amplitude, the growth rate α drops until a limit cycle condition is reached α = 0 for a pressure amplitude | ˆ p | = 2 (βζ ) / (3 κ ) , which is

Fig. 25 .

 25 Fig. 25. (a) Spinning mode calculations. Growth rate α of the two degenerate modes ψ 1 and ψ 2 as a function of the pressure oscillation level | ˆ p | . The limit cycle (black rectangular symbol) is reached for | ˆ p | = Ā = B = 2 (βζ ) / 3 κ. (b) Standing mode calculations. Growth rate α of mode ψ 1 plotted as a function of the maximum pressure oscillation level | ˆ p | . The limit cycle (black rectangular symbol) is reached for | ˆ p | = A = 4 / 3 (βζ ) /κ. Due to symmetry, results corresponding to the distribution ψ 2 are the same. β= 0.15, ζ = 0.05 and κ= 0.2 for both spinning and standing calculations.

  | ˆ p | ( θ ) . As a consequence, modes ψ 1 and ψ 2 are not degenerate. In the calculation, for each pressure level | ˆ p | j , the eigenmode ψ 1 is used to distribute | ˆ p | over θ ( B= 0 and A = 0). Due to symmetry, results with the distribution ψ 2 are the same. Figure 25 b shows results of the simulation for the growth rate α as a function of the maximum pressure oscillation level | ˆ p | max , which is varied from 0 to 1. Analytical results for the growth rate α = 1 / 2(β -ζ ) at vanishing perturbation amplitudes are again captured by the simulation when | ˆ p | max → 0. The growth rate then decays until a limit cycle is reached α= 0 when the maximum pressure oscillation level reaches | ˆ p | max = 4 / 3 (βζ ) /κ, a value again in agreement with

J 11 =

 11 λ 1 = ∂ ˙ A ∂A A= A , B=0 and J 22 = λ 2 = ∂ ˙ B ∂B A= A , B=0 . (17)In these expressions, the dot operator indicates the time derivative. The evaluation of the Jacobian coefficient J 11 starts from the calculation of the time derivative of the mode amplitude A . This can easily be done by recalling the definition of the growth rate ˙ A = αA . This product is calculated for each amplitude using numerical values of the growth rate α and amplitude A determined by the code. Computations of ˙ A (black dashed line) are plotted in

Fig. 26 .

 26 Fig. 26. Numerical evaluation of the Jacobian matrix coefficients for the standing mode. Calculations make use of results for the growth rate α of the mode ψ 1 (black continuous line), for the growth rate α ψ 2 of the ψ 2 mode when the pressure fluc- tuation level | ˆ p | is distributed following the ψ 1 distribution (red continuous line) and ˙ A (black dashed line) as a function of | ˆ p | max . β= 0.15, ζ = 0.05 and κ= 0.2. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

Fig. 26

 26 Fig. 26 together with the growth rate α (black continuous line) as a function of the maximum pressure fluctuation level | ˆ p | max . Following its definition, J 11 corresponds to the slope of the black dashed curve taken at the equilibrium point | ˆ p | max = A . The computed value λ 1 = ζβ coincides with the one predicted analytically [16] . For A → 0, it is also worth noting that the slope of the black dashed curve is equal to the growth rate α = (βζ ) / 2 predicted by a linear stability analysis. The red continuous curve in

Fig. 26

 26 Fig. 26 corresponds to the growth rate α ψ 2 trajectory of the ψ 2 mode when the pressure fluctuation level | ˆ p | is distributed following the ψ 1 distribution. When | ˆ p | max → 0, the black and red trajectories converge to the same value α ψ 2 = α = (βζ ) / 2 and the two modes are again found to be degenerate. Increasing the oscillation level | ˆ p | max , the black and red continuous trajectories diverge. For each point of the red trajectory B= 0, indicating that this mode has no influence on the standing equilibrium point of

Table 2

 2 Synthesis of the stability and trajectory analyses.

		Operating point A	Operating point B
	Stability of spinning mode based	Stable when criterion is applied at | u / u | = 0.61 and	Stable when criterion is applied at | u / u | = 0.59 and
	on Eq. (8)	f = 473 Hz	f = 466 Hz
	Stability of standing mode based	Unstable when criterion is applied at | u / u | = 0.86 and	Stable when criterion is applied at | u / u | = 0.86 and
	on Eq. (10)	f = 473 Hz	f = 478 Hz
	Prediction of stability analysis	Spinning mode	Spinning or standing mode
	Trajectory analysis	Spinning mode features a growth rate that is slightly	Standing mode growth rate is larger than that of the
		greater than that of the standing mode over a limited	spinning mode for all oscillation amplitudes
		range of amplitudes