## Abstract

This theoretical study investigates spinning and standing modes in azimuthally symmetric annular combustion chambers. Both modes are observed in experiments and simulations, and an existing model predicts that spinning modes are the only stable state of the system. We extend this model to take into account the effect that the acoustic azimuthal velocity, *u*, has on the flames, and propose a phenomenological model based on experiments performed on transversely forced flames. This model contains a parameter, *δ*, that quantifies the influence that the transversal excitation has on the fluctuating heat release. For small values of *δ*, spinning modes are the only stable state of the system. In an intermediate range of *δ*, both spinning and standing modes are stable states. For large values of *δ*, standing modes are the only stable state. This study shows that a flame’s response to azimuthal velocity fluctuations plays an important role in determining the type of thermoacoustic oscillations found in annular combustors.

## 1. Introduction

Combustion systems such as aeroplane engines and rocket engines are often susceptible to large amplitude self-sustained pressure oscillations, called thermoacoustic oscillations. These lead to excessive noise and sometimes to structural damage [1]. Self-sustained oscillations occur when the phase difference between pressure fluctuations and heat release fluctuations is less than one-quarter cycle, as described by Rayleigh [2]. These heat release fluctuations arise from the flames’ response to incident velocity or pressure perturbations. The flames’ response therefore plays a crucial role in determining a system’s thermoacoustic behaviour. In this paper, we use a low-order model of a thermoacoustic system to determine the implications of a flame response that has been observed experimentally in annular combustors.

Annular combustion chambers are commonly used in aircraft gas turbines because they fit efficiently between the axial compressor and the turbine. Their circumference is much longer than their length and width, so thermoacoustic oscillations tend to develop in the azimuthal direction. If they travel in a clockwise or anticlockwise direction, with the pressure and velocity nodes travelling at the speed of sound, they are called spinning modes. If the nodes are fixed in space and the wave modulates its amplitude without travelling, then they are called standing modes. Both types of mode are found in large Eddy simulations (LES), experiments and real engines. See, for example, Noiray & Schuermans [3] for spinning modes and Worth & Dawson [4], Wolf *et al.* [5] for both spinning and standing modes.

Schuermans and co-workers [6] study an annular combustor as a network of acoustic elements, using a state space representation. Their linear stability analysis predicts that standing modes are linearly unstable. In time, however, these develop into a spinning mode, which they show is the only stable limit-cycle of the system. They show that this behaviour is also seen for a thermoacoustic model containing a one-dimensional wave equation and a nonlinear saturating pressure-dependent heat release. This model is similar to that which will be used in this paper, in equation (2.11).

Noiray *et al.* [7] consider the effect of a non-uniform heat release in the azimuthal direction. If the acoustic mode has azimuthal dependence of the form , then they show that a non-uniform perturbation of heat release of the form is particularly influential. The amplitude of this non-uniformity is labelled *C*_{2n}. For *C*_{2n}=0, their analysis predicts that only spinning modes are stable. For larger values of *C*_{2n}, a sum of standing and spinning modes can be stable. Above a critical value of *C*_{2n}, only standing modes are stable.

This does not explain, however, why standing modes are the preferred state of the system in some rotationally symmetric configurations with *C*_{2n}=0, as found in Worth & Dawson [4] and Wolf *et al.* [5]. Combustors are very noisy environments, and one explanation could be that noise causes the thermoacoustic oscillations to switch between different modes [8,9]. Noiray & Schuermans [3,10] discuss the effect of noise on the system presented by Noiray *et al.* [7], for a symmetric configuration (by setting *C*_{2n}=0). The only deterministic, stable states of the system remain the two spinning modes, as it is when there is no noise. However, noise can make the system jump between the two modes, and, when it does so, the system passes through the vicinity of a standing mode. Theoretical results and experimental data agree in presenting a probability density function of the state of the system with two clear peaks on the two spinning modes. However, this is not consistent with Wolf *et al.* [5] and with certain configurations of Worth & Dawson [4], where the system has a statistical preference for standing modes.

In summary, current thermoacoustic models cannot explain why standing modes in symmetric annular chambers should be a preferred state of the system, despite experimental evidence that they sometimes are. In this paper, we extend the work by Noiray *et al.* [7] to include the influence of transversal flame excitation and show that a phenomenological model that includes transverse excitation can exhibit stable standing modes as well as stable spinning modes.

The paper is organized as follows: in §2, we present a concise derivation of the one-dimensional equation (2.11) governing the problem. In §3, we discuss the model for the heat release, taking into account the transversal forcing in §4. In §5, we discuss how to simulate equation (2.11) numerically and present some introductory results. In §6, we reduce the problem to a system of coupled oscillators, and in §7, we discuss its stability. We then apply the method of averaging and study the resulting phase space, providing a graphical description of the system in three dimensions.

## 2. Description of the problem

The geometry under investigation is a thin annular combustion chamber. We study this problem in cylindrical coordinates and time: . *z** is the longitudinal length of the combustor, and *R*_{1} and *R*_{2} are the radii of the inner and outer surfaces. For the sake of simplicity, we consider *R*_{1} and *R*_{2} to be uniform in *z*. We consider mean and fluctuating variables, neglecting the influence of viscosity on the flow field. The momentum and pressure equations for the fluctuating variables, considering only first-order acoustics [11], are
2.1and
2.2In the above expressions, the prime indicates fluctuating components, and are the mean density of the gas mixture and the mean pressure, *γ*=*c*_{p}/*c*_{v} is the ratio of specific heat capacities, and *q*′ is the fluctuating heat released by the combustion. Spatial averaging has already been applied, and the effects of a non-uniform speed of sound have been neglected [11]. This is common to many other studies modelling this problem [7,12]. Terms of order are neglected by assuming a low Mach number flow.^{1} We study the problem in polar coordinates and drop the dependence on the radial coordinate because this is weak, even if the gap between the two cylinders is non-negligible [13]:
2.3and
2.4Here, *u*′ is the component of the velocity in the azimuthal direction and *R*≈(*R*_{1}+*R*_{2})/2. We proceed by non-dimensionalizing the equations, picking a new time scale:
2.5where is the spatially averaged speed of sound. We introduce the non-dimensional variables *p**,*q** and *u**, defined as
2.6
2.7
and
2.8The non-dimensional system of equations, dropping the asterisks and the prime on the time variable, becomes
2.9and
2.10In this expression, we have included a damping term, with a coefficient *α*>0. Equations (2.9) and (2.10) are equivalent to the wave equation, with ∂*q*/∂*t* as a source term:
2.11This model has been the common starting point of Noiray *et al.* [7] and Noiray & Schuermans [3,10]. It is the one-dimensional counterpart of the model studied by Sensiau *et al.* [12], where the whole three-dimensional field is considered. This simplification is not appropriate for the study of radial and longitudinal instabilities, but allows an analytical treatment of azimuthal instabilities.

## 3. Heat release constitutive equation

The heat release fluctuations, *q*, are often assumed to depend either on velocity, pressure, or both, with the inclusion of one or more time delays. However, we start from the analysis of Noiray *et al.* [7], where *q* is a function of *p* only:
3.1and
3.2The first term depicts linear growth governed by *β* for small fluctuating pressures. The second term in (3.2) is a nonlinear cubic saturation, governed by a coefficient *κ*. For any given *κ*>0, the study of equation (2.11) in terms of a new pressure variable leads to a new problem independent of *κ*. It follows that the coefficient *κ* induces simply a rescaling of the problem, and will be set to 1 in the following analysis. By assuming *q*=*f*(*p*) to be an odd function of *p*, (3.2) is a fourth-order-accurate Taylor expansion of *f*. The analysis carried out in this paper could easily be extended to higher-order terms.

The theory developed in Noiray *et al.* [7] based on (3.1) does not predict stable standing modes for symmetric configurations, which are observed as preferred state of the system in Worth & Dawson [4] and Wolf *et al.* [5]. The universal validity of (3.1) is then called into question, particularly the assumptions of (3.1), which are (i) the absence of a time delay in *p* and (ii) the independence of *q* on anything else except *p*.

Regarding the first point, *q* has been found to be reasonably in phase with *p* in an LES simulation of a specific, symmetric rig [5]. In that rig, both standing and spinning modes are observed, suggesting that standing modes are possible in the absence of a time delay. Because the aim of this paper is to explain how standing modes are possible in symmetric systems, we do not consider a time delay and assume that *p* and *q* are in phase, leaving this investigation for further research. Some results based on linear stability, applied to an *n*−*τ* model, discussing the importance of the time delay can be found in Sensiau *et al.* [12].

Regarding the second point, one possibility is to assume that the heat release depends on the azimuthal coordinate, because combustion happens mainly near the injectors. This can be performed by introducing a shape function *ψ*(*θ*), which is large near the flames and small far from the flames:
3.3and
3.4where *M* is the number of injectors.

The expression (3.4) is one of the many possibilities for such a shape function. The adoption of the constitutive equation (3.3) instead of (3.1) in the analysis we will develop next does not lead to any qualitative differences in the stability analysis of the standing and spinning modes, and will therefore not be discussed further. This means that a spatially accurate description of the heat release does not explain standing modes in symmetric annular chambers, and therefore suggests that it is sufficient to consider spatially averaged models for the heat release *q*. We stress that this statement does not apply to non-symmetric annular combustors, in which the shape function *ψ* is no longer *M*-periodic in the azimuthal direction, and in the very specific case where *M*=2*n* in (3.4), where *n* is the unstable acoustic mode. This latter case is covered by the stability analysis of Noiray *et al.* [7], and shows that mixed modes and standing modes are possible. In this case, the shape function has peaks at the pressure antinodes, and troughs at the pressure nodes of the standing mode.

One other possibility is to assume, in addition to the dependence of *p*, a dependence of *q* on the azimuthal velocity *u*, which excites transversally the flames. This possibility is investigated in this paper.

## 4. Model of transversal forcing

The effect of transverse excitation on swirling premixed flames is a current topic of research. Hauser *et al*. [14] report that an asymmetric perturbation of higher hydroxil radical intensity is generated by the transversal velocity. This asymmetric region of stronger combustion spirals around the injector at the forcing frequency. This asymmetry also persists in addition to longitudinal forcing, suggesting that the two phenomena are superposable.

We report in figure 1 two phase-averaged chemiluminescence images of flames at a pressure antinode (velocity node) and at a pressure node (velocity antinode), kindly provided by James Dawson and Nicholas Worth.

In figure 1, at pressure antinodes, there is no transverse velocity excitation. Circles of positive/negative heat release are shed from the injector and propagate outwards. The fluctuating heat release is found to be approximately axisymmetric around the injector: at every instant in time of a limit-cycle, the phase of the perturbation is approximately axisymmetric. This is consistent with O’Connor & Lieuwen [16], where the vorticity disturbance is symmetric around the injector at pressure antinodes. At velocity antinodes, the symmetry of the perturbation breaks: the heat release is found to be approximately in anti-phase on the two sides of the flame, in the direction of the transverse velocity. In O’Connor & Lieuwen [16], the same break of symmetry happens for vorticity disturbances, which are asymmetric at velocity antinodes. This means that the spatially averaged heat release fluctuation of an injector is smaller at velocity antinodes. This happens because, in the averaging, the zones in anti-phase cancel out in figure 1*b*. Because the flame diameter is small compared with the wavelength of the unstable acoustic mode, we can assume that the flame is acoustically compact, and spatially average the heat release on each burner.^{2} Based on this observation, we assume that the fluctuating heat release of an injector is smaller if a transverse excitation is present, introducing a dependence on the velocity *u*:
4.1In expression (4.1), *f* is the same function introduced in (3.2), and all the previous considerations apply to it. The function *μ* is the degree of symmetry of *q* around a burner, and must be unity for zero transverse excitation, and smaller than 1 for *u*≠0, in the range of velocities investigated:
4.2We study two possibilities for *μ*:
4.3and
4.4where *δ* is a positive coefficient that expresses how strongly the transverse forcing influences the heat release. Both models must respect (4.2) at every instant of time. The structure (4.1) and the choice of *μ* in (4.3) and (4.4) are not intended to be an accurate representation of the flame’s behaviour, because they are based on qualitative observations. The model is phenomenological, and the aim is to discuss the effect of transverse forcing on the stability of standing and spinning modes on a qualitative level.

## 5. Numerical set-up

Because it is not straightforward to evaluate the time derivative of *q* in equation (2.11), we opt to study the system of equations (2.9) and (2.10), where such a derivative is not required. We project the equations into Fourier space, obtaining a system of ordinary differential equations. The generic *n*th complex Fourier mode is governed by the equations
5.1Because *q* is a nonlinear function of *u* and *p*, at each timestep the two functions *u*(*θ*) and *p*(*θ*) are evaluated from the Fourier coefficients {*u*_{n}} and {*p*_{n}}, and then *q*(*θ*) is calculated as *f*(*p*(*θ*))*μ*(*u*(*θ*)). Finally, the {*q*_{n}} coefficients are evaluated as a Fourier transform of *q*(*θ*). The system (5.1) can then be numerically integrated with a numerical scheme.

The damping of this problem has to be adjusted to avoid excessive growth of higher-order harmonics.

Specifically, we consider only the dissipation owing to the boundary layers, which scales as the square root of the frequency [11,17]. We take this into account fixing in (5.1).

Two examples of two simulations showing a spinning and a standing mode are reported in figure 2, truncating the number of Fourier modes to *N*_{f}=161. The two pictures do not imply that the two modes are stable, and only time marching for a long time allows us to check this at this stage.^{3} We present here both cases A and B only to show what the two different *μ* functions defined in (4.3) and (4.4) look like. Both cases present spinning modes for small values of *δ* and standing modes for large values of *δ*.

The existence of a standing mode at one value of *δ* and a spinning mode at another is a key result of this paper. In the following sections, we conduct a stability analysis of these modes to confirm that they are indeed both stable limit cycles of the nonlinear governing equations.

## 6. Reduction to a system of coupled oscillators

In this section, we carry out spatial averaging in the azimuthal direction, in the same way carried out in Noiray *et al.* [7]. When annular combustors are subject to azimuthal instabilities, there is usually only one strong Fourier component, which corresponds to the *n*th lowest acoustic mode of the chamber. This is apparent from the power spectral density of the Fourier transform of pressure signals from experiments Worth & Dawson [15], and is also observed in the numerical solutions of (5.1). We truncate the modal expansion and consider only the *n*th mode:
6.1and
6.2where the second expression was obtained substituting (6.1) into (2.9). We now apply spatial averaging [11] to this system: we substitute (6.1) and (6.2) into (2.10), multiply the expression by , and then average over 2*π* in the azimuthal coordinate, obtaining (6.3):
6.3and
6.4Here, (6.4) was obtained similarly by multiplying by . Note that these expressions are exact, and the assumption that higher-order modes are negligible is applied assuming that the two source terms *F*_{i} on the RHS depend only on the Fourier modes *η*_{1} and *η*_{2}. They are
6.5and
6.6We can study the system in the new time scale *t*′=*nt*, and obtain
6.7and
6.8where *α*↦*α*/*n*, and the expressions of *f*_{i}≡*F*_{i}/*n*^{2} and how to evaluate them are reported in appendix A. This is a system of coupled oscillators, which can be numerically integrated in time in a four-dimensional phase space, as opposed to the phase space with *N*_{f} dimensions introduced in (5.1).

### (a) Amplitudes and phase representation

Instead of studying the system in terms of displacements, *η*_{i}, and velocities, *η*_{i}′, it is more useful to study it in terms of amplitudes and phases,
6.9
6.10
6.11
and
6.12In these expressions, the frequency *ω* can, in principle, be perturbed by the parameters of the problem from the natural frequency of the two oscillators, which is 1. Equations (6.7) and (6.8) are symmetric in *η*_{1},*η*_{2}. It is useful to introduce the phase difference between the two oscillators, *ϕ*(*t*)≡*φ*_{1}(*t*)−*φ*_{2}(*t*). If *ϕ* settles to ±*π*/2 and *A*=*B*, then the substitution of (6.10) and (6.12) into (6.2) shows that the pressure distribution corresponds to a spinning mode in the counterclockwise/clockwise direction, respectively:
6.13This solution spins in the azimuthal direction as in figure 2*a*. On the other hand, if *ϕ* settles to *π* or 0, then there is a standing mode, for any value of *A*,*B*:
6.14The pressure nodes can be found by studying the zeros of the *θ*-term in (6.14). They are fixed in space, as shown in figure 2*b*. It is convenient to examine the two cases in the (*η*_{1},*η*_{2}) plane as a function of time. With reference to equations (6.9)–(6.12), the two modes give rise to limit-cycles which are either circles or lines. The situation is summarized in table 1. Figure 3 shows two simulations of trajectories in the (*η*_{1},*η*_{2}) plane for case B. This is similar to fig. 11 in Paschereit *et al.* [6]. It is worth noting that the complex number *C*(*t*)=2(*η*_{1}(*t*)+*i* *η*_{2}(*t*)) is the indicator proposed by Poinsot *et al.* [8] to study the nature of these modes. The two cases have different values of *δ*, and lead to either spinning or standing limit-cycles. The simulations have been started with nearly the same initial condition for *ϕ*, and from two random values for *A*,*B*.

The main objective of this section is to reduce the original partial differential equation to a system of coupled oscillators, and to present a simpler way to look at standing and spinning modes in terms of amplitudes and phase in the (*η*_{1},*η*_{2}) plane. In §7, we will perform a stability analysis of these modes, and in §9, we will present a phase space realization of the system in terms of the two amplitudes and of the phase here introduced.

## 7. Stability of the coupled oscillator system

We first report some results from the linear analysis of the fixed point *p*(*t*,*θ*)=*u*(*t*,*θ*)=0 of the system (6.7) and (6.8). This fixed point is stable for *β*<*α*. A double Hopf bifurcation occurs at *β*=*α*, where two complex eigenvalues cross the imaginary axis at the same time. Similarly, Sensiau *et al.* [12] perform a linear stability analysis of an azimuthally symmetric chamber and find two linearly unstable spinning modes with exactly the same growth rate. They conclude that, in perfectly symmetric systems, the sum of the two identical spinning modes would lead to a stable standing mode. This is not the case, however, as shown by the fact that, for *δ*=0, and *β*>*α*, this system converges to a stable spinning mode, in accordance with Paschereit *et al.* [6].

We proceed by analysing the case of the oscillating system, fixing *α*=0.08 and *β*=0.10, and focusing on case A. Figure 3 shows that two different values of *δ* lead to two different limit cycles: a spinning mode and a standing mode. We now study the system over a range of *δ*. To do this, we numerically integrate the system until it converges to a limit cycle, and then track the limit cycle as we vary *δ* using MatCont, a numerical continuation package [18]. Figure 4 shows the stability of the spinning and standing modes.

For *δ*=0, the spinning mode is stable, because all its Floquet multipliers are smaller than 1 in figure 4*a*. At *δ*_{c2}≈1.027, the modulus of two Floquet multipliers crosses 1, which corresponds to a subcritical^{4} Neimark–Sacker bifurcation at which the spinning mode becomes unstable. The argument of these two Floquet multipliers, shown in figure 4*c*, is small. Note that, for *δ*>*δ*_{c2}, there is one multiplier with modulus smaller than 1 and two multipliers with modulus greater than 1. This means that the system is attracting from an invariant manifold^{5} with dimension 1, and repelling to another invariant manifold with dimension 2. This is consistent with figure 3*b*, where the point is first attracted to the spinning mode (circular line) before being repelled towards a standing mode (straight line).

For the standing mode, for every value of *δ*, two multipliers are exactly equal to 1. One of these is due to the fact that the system is at a limit-cycle and any movement in the direction of the limit-cycle remains on the limit-cycle (the spinning mode has one too, under the horizontal black line in figure 4*a*). The other is due to the fact that the nodes of the standing mode can rotate arbitrarily around the annulus, i.e. the black line in figure 3*b* can take any angle with the axes. A fold bifurcation occurs at *δ*_{c1}=0.949, making the mode stable for *δ*>*δ*_{c1}. Note that, for *δ*<*δ*_{c1}, there is one multiplier with modulus smaller than 1 and one multiplier larger than 1. This means that the system is attracting from a one-dimensional invariant manifold and repelling from another one-dimensional invariant manifold. This can be seen in figure 3*a*, where for a while the solution lingers as a standing mode (straight line) before being repelled towards the spinning mode (circle line).

The angular frequency of the limit cycles is not changed by the nonlinearities of the problem: in the range of parameters investigated, the period of oscillations was found to be constant and equal to 2*π*. In summary

— for

*δ*<*δ*_{c1}≈0.949, only the spinning mode is stable;— for

*δ*>*δ*_{c2}≈1.027, only the standing mode is stable;— for

*δ*_{c1}<*δ*<*δ*_{c2}, the system is multistable, with both standing and spinning modes stable.

Moreover, for *δ*<*δ*_{c1} and *δ*>*δ*_{c2}, the unstable mode attracts the solution on a one-dimensional invariant manifold, before repelling it towards the stable mode.

We checked that these stability results, obtained for the system (6.7) and (6.8), apply also to the original system (2.11) by performing numerical simulations of (5.1) for different values of *δ*. For each value of *δ*, we started the simulation with both standing and spinning modes as initial conditions, and evaluated their stability by time-marching. The same qualitative picture was found, with the two critical values of *δ* confined in these intervals: 0.9<*δ*_{c1}<1.0 and 1.1<*δ*_{c2}<1.2, in good agreement with the values just presented. This shows that the reduction to a system of coupled oscillators by considering only the fundamental unstable harmonic, as presented in §6, is a powerful tool to study the stability of the original wave equation (2.11), at least for the values of *α*,*β* investigated here.

We do not report here the analysis for case B, because the overall behaviour is the same as that of case A.

## 8. Slow flow

In this section, we apply the method of averaging to the system of coupled oscillators (6.7) and (6.8) for case B. We will obtain a new system of differential equations in terms of the amplitudes of oscillation *A*,*B* and of the phase difference *ϕ*, introduced from equation (6.9) onwards. This will reduce the dimensions of the problem from 4 to 3, allowing us to visualize the complete dynamics of the problem. The method of averaging [20] gives the following formulation of the slow flow:
8.1where and , and the averaging operator of a generic function *h* is introduced as
8.2Note that, while in the definitions (6.9)–(6.12), the amplitudes and the phases are functions of time, they are constants in the RHS of (8.2). We fix *ω*=1, consistent with the period being 2*π* as reported earlier. Some details on how to tackle the four integrals can be found in appendix B, together with the full equations of the system (8.1).

From now on we fix, as previously, *κ*=1,*α*=0.08,*β*=0.1. The two critical values of *δ* for standing and spinning modes are respectively *δ*_{c1}≈6.2076 and *δ*_{c2}≈6.2165. We can then visualize this phase space in terms of *A*,*B* and *ϕ*, as a function of *δ*. The amplitudes *A*,*B* are non-negative numbers, and *ϕ*∈[0,2*π*]. Because the phase space is symmetric with respect to the planes *ϕ*=*kπ*/2 with *k*=0,1,2, we restrict the visualization to *ϕ*∈[*π*/2,*π*]. The system is also symmetric with respect to the plane defined by *A*=*B*.

Because of the difficulty of drawing a three-dimensional phase space, we report the flow on a few invariant manifolds. These completely describe the stability of the problem.^{6} Figure 5 shows two convenient slices of the same phase space for *δ*=3. In figure 5, every shaded surface is an invariant manifold, and all invariant manifolds are reported, with the exception of the two planes *A*=0 and *B*=0, and the plane *A*=*B*, which is reported in figure 6. In figure 5, only the spinning mode is stable, because *δ*<*δ*_{c1}. We then fix *δ*=12>*δ*_{c2} and present the same slices of the phase space in figure 7, in which only the standing mode is stable.

In this representation, the addition of a non-zero asymmetry parameter, *C*_{2n}, as proposed in Noiray *et al.* [7], shifts the red point of figure 5 towards one of the *A*,*B* axes, maintaining it on the same plane *ϕ*±*π*/2. Doing so, the system exhibits a superposition of standing and spinning modes. As discussed in Noiray *et al.* [7], above a certain threshold, the red point hits and gets stuck on one of the *A*,*B* axes, becoming a pure standing mode. The current analysis shows that, with the addition of transverse forcing introduced in (4.1), with *μ*(*u*) from case B as defined in (4.4), the standing mode becomes stable in a different way, without passing through a superposition of standing and spinning modes. It is worth noting that this analysis could easily be extended to include the parameter *C*_{2n} in order to discuss the stability of transverse forcing in non-symmetric annular chambers.

## 9. Conclusions

This study improves the current understanding of standing and spinning modes in symmetric annular combustion chambers, which is the subject of current research [5,6,10,12,15]. The starting point of this study is the model proposed in Noiray *et al.* [7]: the fluctuating heat release *q* is assumed to grow linearly and saturate nonlinearly as the pressure increases, as *q*=*f*(*p*)=*β*−*κp*^{3}. In our analysis, we add an extra dependence, which reflects experimental observations [14]: the fluctuating heat release fluctuates axisymmetrically at velocity nodes (pressure antinodes), whereas it fluctuates from side to side at velocity antinodes (pressure nodes), as in figure 1. When integrated over a sector of the chamber, the *q* fluctuations are larger in the first case. We then assume that *q*=*f*(*p*)*μ*(*u*), and we study two ways in which *q* can depend on *u*. We consider a case A with *μ*(*u*)=1−*δ*|*u*|, and a case B with *μ*(*u*)=1−*δu*^{2}.

For both cases, we find that (i) for small *δ*, only spinning modes are stable; (ii) for intermediate *δ*, both standing and spinning modes are stable, and the system is multistable; (iii) for large *δ*, only standing modes are stable. We show that this standing mode is fundamentally different from the one found in Noiray *et al.* [7] in non-symmetric chambers, and it affects the phase space in a different way, as described in §8 for case B.

Another result is that, when the system has only one stable limit-cycle, the other unstable limit-cycle is not a repellor: it attracts the solution on one invariant manifold, and repels it on another. Figure 12 in Paschereit *et al.* [6] suggests that the same property also applies to their system. If this property holds in industrial combustors, then noise could randomly shift the point in the phase space also to the attracting manifold of the unstable mode, and the system could linger for longer close to the unstable mode before decaying to the stable one. We give an example of this transient behaviour in figure 3, and we comment on it based on the stability results in §7.

The analysis can be extended by adding higher-order terms in *u* and *p* to the model for *q* in order to match the model to results from experiments on a single injector. The analysis can be extended to complex geometries, as long as the flame is acoustically compact.

We mention that the oscillators formulation (6.7) and (6.8), and the slow flow equations (8.1) can both be used as a physically based model of the combustion process for purposes of control.

This study suggests that transversal forcing plays an important role in annular combustion instabilities, and should be taken into account to accurately predict instabilities in annular configurations. The experimental characterization of a single injector to longitudinal forcing seems to not be sufficient to predict the final state of the combustor.

## Funding statement

This work was supported by the European Research Council through project ALORS N.259620.

## Acknowledgements

We are very grateful to Nicholas Worth and James Dawson for making their data available at a previous stage of this work, for offering figure 1, and for helpful discussions on their experiment.

## Appendix A

**(a) Spatial averaging**

Here, we include the terms on the RHS of (6.7) and (6.8), evaluated with the help of a computer algebra system. For case B, they are
A1and
A2For case A, the evaluation of the integrals (6.5) and (6.6) is difficult. It is necessary to split the integrals into two domains, in which the argument of the absolute value is either positive or negative, and then put together the results afterwards. The boundaries of these domains depend on the argument *φ* of the complex number *η*_{1}+*iη*_{2}. The final result is therefore in terms of *φ*, and it is too long to be written here. It depends on cosine and sine functions of *ϕ* and its multiples. It is possible to eliminate the dependence on *φ* by trigonometrically expanding the terms, and substituting these relations for the sine and cosine of *φ*:
A3and
A4Both cases A and B contain forcing terms that are the same as the ones obtained in Noiray *et al.* [7] for *δ*=0.

**(b) Method of averaging**

Here, we show how to evaluate the terms of (8.1) for case B. It is possible to prove that
A5where we did not indicate the explicit dependence of *f*_{1} on *A*,*B*,*θ* as in the definition (8.2) of the averaging operator for conciseness, and we set *ω*=1. Then, from *F*_{1}, we can evaluate the two terms on the left as real and imaginary parts. We can set *z*=e^{it}, and operate these substitution in *F*_{1}:
A6

We can then change the line integral in *F*_{1} to a contour integral on the unit circle of the complex plane:
A7It can be shown that this function presents a single pole at the origin with no branch cuts, so that
A8Note that the expression for *f*_{1} is long, and the evaluation of the residue requires the use of a computer algebra system. The system of equations of the slow flow is
A9
A10
and
A11

## Footnotes

↵1 The azimuthal component of the mean field is induced only by the injectors’ swirlers, and is negligible.

↵2 See also the comment on the use of spatially averaged models for

*q*in §3.↵3 The two modes in the two cases will be proved to be stable later with rigour.

↵4 Based on the first Lyapunov exponent, which is positive.

↵5 For the purposes of this paper, an invariant manifold can be thought of as a particular surface in the phase space such that all points on it are either attracted to or repelled from the same limit-cycle or fixed point. Refer to [19], for a rigorous definition.

↵6 The flow perpendicular to an invariant manifold is zero.

- Received April 12, 2013.
- Accepted June 13, 2013.

- © 2013 The Author(s) Published by the Royal Society. All rights reserved.