Linear Stability Analysis of Transversal Thermoacoustic Modes in Reheat Combustors

Simon M. Heinzmann 1 Institute of Energy Systems and Fluid-Engineering, ZHAW Zurich University of Applied Sciences, Winterthur 8401, Switzerland; Department of Mechanical and Process Engineering, ETH Z€ urich, Z€ urich 8092, Switzerland e-mail: [email protected] Harish S. Gopalakrishnan Institute of Energy Systems and Fluid-Engineering, ZHAW Zurich University of Applied Sciences, Winterthur 8401, Switzerland Birute Wood Ansaldo Energia Switzerland, Baden 5400, Switzerland Mirko R. Bothien Institute of Energy Systems and Fluid-Engineering, ZHAW Zurich University of Applied Sciences, Winterthur 8401, Switzerland Linear Stability Analysis of Transversal Thermoacoustic Modes in Reheat Combustors Gas turbines featuring sequential combustion architectures are ideally suited to deliver CO2-free, fuel-flexible, and on-demand power to the grid. A sequential combustor is comprised of two axially staged flames. Typically, the first stage is a swirl-stabilized propagating flame. The hot products of the first stage are diluted with air before additional fuel is injected in the second stage. The resulting vitiated mixture leads to auto-ignition in the second combustor. The sequential architecture can be leveraged to vary the fuel split between the combustion chambers to account for different reactivities of alternative fuels. Thermoacoustically stable combustors are critical to ensure low emissions, high reliability, and mechanical integrity. Under specific off-design conditions, the auto-ignition-stabilized second stage can be subject to transversal instabilities. This article presents a finite-elementcoupled method to model the thermoacoustic behavior of the Ansaldo Energia H-class GT36 reheat combustion stage in a very cost-effective manner. To do so, the auto-ignition flame is divided into multiple parts for which the flame transfer function (FTF) methodology commonly used for planar waves is applied. The stability of transverse eigenmodes is assessed and validated against experiments performed under engine conditions. We show that the framework can correctly predict the stability of high-frequency combustor modes. Consequently, it can be used to get reliable stability estimates of high-frequency thermoacoustic modes in industrial reheat combustion chambers and to ensure that such instabilities are avoided. [DOI: 10.1115/1.4069453] Keywords: thermoacoustics, transverse modes, high-frequency, reheat combustion, stability analysis 1 Introduction To counteract global warming, a rapid reduction of green-house gas emissions is needed. The electricity sector is foreseen to lead the way in terms of net zero emissions by 2050 [1]. In order to achieve this, the future energy production will strongly rely on renewable energy systems such as solar-, wind-, or hydropower. However, these sustainable energy solutions are subject to fluctuating electricity production on a broad time-scale [1]. To balance the intermittent nature of energy production and consumption, excess power produced from renewable energy systems can be chemically stored in renewable fuels such as hydrogen, ammonia, or synthetically produced conventional fuels. Subsequently, these fuels can be fired in gas turbines to produce environmentally friendly ondemand power, making them an ideal asset in the future energy landscape [2]. To fire a variety of renewable fuels in gas turbines, the constant pressure sequential combustion architecture, as implemented in the Ansaldo Energia GT36, is a promising solution [2–4]. A sequential combustor consists of two combustion chambers arranged in series (Fig. 1). The multiburner first stage consists of multiple propagationstabilized burners [4]. Downstream of the first stage, dilution air is admixed in the dilution air mixer. Further downstream, fuel is admixed in the sequential stage center body burner (CBB), and due to the high temperature of the resulting mixture, an auto-ignitionstabilized flame, also called reheat flame, results [2,3,5]. To enable the combustion of fuels with different reactivities, the fuel and air split can be varied between the two combustion stages. This allows to tune the combustion process in each stage separately to account for the different reactivities of gas turbine fuels, thereby expanding the operational flexibility of the engine [5–7]. Thermoacoustically stable combustion stages are important to ensure low emissions, high reliability, and mechanical integrity. Multiple methods exist to perform thermoacoustic analysis of combustion chambers. High-fidelity methods like direct numerical simulations and large eddy simulations provide results with high accuracy. However, simulating a full-scale industry combustor for all relevant operating conditions remains computationally far too expensive. Computationally more affordable modeling techniques on the other hand have shown to also reproduce thermoacoustic behaviors accurately. These methods are typically based on finiteelement (FE) modeling techniques or network-modeling approaches [8–11]. Nicoud et al. [8] showed that the finite-element-method (FEM) can be used to classify the stability of thermoacoustic longitudinal eigenmodes in two-dimensional ducts using a Turbo Expo: Turbomachinery Technical Conference & Exposition (GT2025), June 16–20, 2025. GT2025. 1 Corresponding author. Manuscript received July 4, 2025; final manuscript received July 21, 2025; published online September 24, 2025. Editor: Jerzy T. Sawicki. Journal of Engineering for Gas Turbines and Power DECEMBER 2025, Vol. 147 / 121019-1 Copyright V C 2025 by ASME; reuse license CC-BY 4.0 Downloaded from http://asmedigitalcollection.asme.org/gasturbinespower/article-pdf/147/12/121019/7530914/gtp-25-1285.pdf by guest on 06 October 2025 Helmholtz solver. Laera et al. [9] accurately calculated the growth rates for spinning and standing azimuthal modes of the laboratoryscale MICCA combustion chamber using a 3D Helmholtz solver. Bothien et al. [10] and Schuermans et al. [11] demonstrated that network models can predict the stability of real engine combustors for which the compact flame assumption is satisfied. Hummel et al. [12] coupled a FEM simulation to a network model to predict the stability of transverse modes in a combustor with a propagationstabilized flame. To the best of our knowledge, prediction of growth rates of transverse modes in auto-ignition-stabilized combustors has not been reported in the literature until this point. The principal goal of this study is to develop a framework that is capable of predicting the linear stability of transverse thermoacoustic modes in a realistic gas turbine combustor. To achieve this objective, we employ a FEM-based approach to model the combustor acoustics along with a recently extended Lagrangian particle-tracking-based framework [13,14] to model the autoignition-stabilized flame dynamics. The framework is validated against experimental data of the Ansaldo Energia GT36 reheat combustor at two operating points (OP). 2 Validation Data From GT36 Combustor In this article, the GT36 sequential combustion chamber is considered (see Fig. 1), and the developed stability prediction framework is aimed to reproduce the stability behavior of this combustor. The sequential liner connecting the CBB with the turbine inlet is simplified in order to be rotationally symmetric. Its radius is adapted accordingly to reproduce the same axial change in cross-sectional area as the real engine geometry (Fig. 1). It is verified that the influence on the investigated acoustic modes is negligible in terms of modeshape and frequency. The thermoacoustic simulations are performed for two distinct operating points at engine condition. OP1 corresponds to a thermoacoustically unstable combustor oscillating intermittently. At OP2, the combustor is stable. The difference in temperature between the two OPs was chosen to be larger than 100 K, which also indicates a large stable operating window. Figure 2shows the probability density function (PDF) of the pressure PðgÞin conjunction with the time signal as well as the frequency spectrum for OP1 (Fig. 2(a)) and OP2 (Fig. 2(b)). The intermittent nature of the pulsations is clearly visible in the time signal of the pressure for OP1. The frequency spectrum shows that the main pulsations occur at a normalized frequency of 1 (all frequency values in this paper are normalized using the value of the main frequency peak of OP1). 3 Framework Development The framework presented in this work was initially developed in two dimensions (2D) and validated against measurements of thermoacoustic oscillations observed in a laboratory-scale reheat combustor [14]. It was shown to be able to predict the growth rate of the dominant unstable mode. Comparison to chemiluminescence imaging showed a qualitatively accurate representation of the spatial heat release rate (HRR) perturbations [14]. Here, we expand the framework to three dimensions (3D) to model transverse eigenmodes of the GT36 CBB. The framework components are: a reactive Reynolds-averaged Navier–Stokes (RANS) simulation to obtain the time-averaged base state for stability analysis, a 1D Lagrangian particle tracking framework to compute the flame response [13,14], and a 3D flame segmentation method to account for the acoustic noncompactness of the flame in the transversal direction. Further, an analytical HRR perturbation formulation in 3D and a 3D Helmholtz solver to compute the eigenvalues are used. 3.1 Numerical Setup and Results of the Three-Dimensional Reynolds-Averaged Navier–Stokes Computational Fluid Dynamics Simulations. In this section, the setup of the computational fluid dynamics (CFD) simulation is discussed, and insights into the results thereof are given. The mean temperature and HRR fields are computed with reactive RANS simulations using the commercial software FLUENT. The combustion is modeled using an in-house reheat combustion model based on a progress variable approach and a tabulated chemical source term. The presumed PDF method is used for the turbulence–chemistry interaction. The model is implemented via the FLUENT user-defined function. The chemical source term is tabulated as a function of fuel and air mixture fractions and normalized progress variable dimensions. The slow chemistry is represented via a composite progress variable, which is defined as the sum of intermediate and product species. Fuel, cooling air mixture fractions and unnormalized progress variable transport equations are solved as standard FLUENT userdefined scalar equations. The first moment (mean mixture fractions) for Zf(fuel), ZCA (cooling air), and Yc(unnormalized progress variable) are solved using Eq. (1) [5], where uiis the velocity component, xiis the spatial coordinate, and Ceff is the effective diffusion coefficient Fig. 1 Ansaldo Energia GT36 constant pressure sequential combustion chamber schematic from Ref. [4] 121019-2 / Vol. 147, DECEMBER 2025 Transactions of the ASME Downloaded from http://asmedigitalcollection.asme.org/gasturbinespower/article-pdf/147/12/121019/7530914/gtp-25-1285.pdf by guest on 06 October 2025 @� q~ Zf @tþ@� q~ui~ Zf @xi − @ @xi Ceff @~ Zf @xi !¼0 @� q~ ZCA @tþ@� q~ui~ ZCA @xi − @ @xi Ceff @~ ZCA @xi !¼0 @� q~ YC @tþ@� q~ui~ YC @xi − @ @xi Ceff @~ YC @xi !¼�qe_ xc (1) The second moments Z00 (mixture fraction variances) are solved using the following equation [5]: @�q~ Z00 f @tþ@�q~ui~ Z00 f @xi − @ @xi Ceff @~ Z00 f @xi !¼2:86lturb @~ Zf @xi !2 − 2�q~ Z00 f sturb @� q~ Z00 CA @tþ@� q~ui~ Z00 CA @xi − @ @xi Ceff @~ Z00 CA @xi !¼2:86lturb @~ ZCA @xi !2 − 2� q~ Z00 CA sturb @�q~ Y00 c @tþ@�q~ui~ Y00 c @xi − @ @xi Ceff @~ Y00 c @xi !¼2:86lturb @~ Yc @xi !2 − 2�q~ Y00 c sturb (2) The turbulence time scale sturb is modeled as the scaled specific rate of dissipation ð0:09xGEKOÞ. Equations for Zfand ZCA are conserved scalar equations without a chemical source term. The transport equation for Yccontains a chemical source term, which is closed with a 2D presumed beta PDF for the fuel and cooling air mixture fractions and 1D presumed beta PDF for the normalized progress variable (C). The normalized progress variable is computed as the ratio of the local progress variable Ycand its corresponding equilibrium value Yðc,eqÞ. The chemical source term is given by the following equation [5]: e_ xC¼ððð1 0 _ xZf,ZCA,C ð ÞPDF Zf,ZCA ð ÞPDF C ð ÞdZfdZCA dC (3) In the 3D numerical analysis, one-twelfth of the burner sector (spanning one full single injection finger) is considered (Fig. 3(a)). The combustor sector transition section volume is matched to the actual value of the engine geometry. The domain is meshed with a tetrahedral mesh with prism boundary layers at the walls. The mesh size is approximately 20 �10 6 cells with refinements applied in the burner mixing section and expected flame location. The mesh resolution is sufficient to resolve the flame brush with multiple cells, see Fig. 3(b). The near-wall regions are meshed y þ <15.0. The steady 3D RANS simulations are performed using FLUENT’s GEKO model (generalized k–xGEKO model) with tuned CMIX and CSEP parameters (“mixing” and “separation” parameters in GEKO model). Fuel is modeled as natural gas with 7% C 2 þ(93.11% CH 4 , 6.51% C 2 H 6 , and 0.38% inerts by volume). The hot gas inlet is assumed to be at equilibrium composition of fuel–air mixture at target hot gas inlet temperature. A circumferentially averaged temperature profile is applied at the combustor inlet. The temperature profile was obtained from first stage combustor and dilution air mixer simulations. No radiation heat transfer is included. A coupled scheme is applied to pressure–velocity coupling with the PRESTO! (PREssure STaggering Options) scheme for pressure. All equations are discretized using a second order upwind scheme (mass, momentum, energy, species, and user scalars). The temperature field on the centerplane shows that OP1 (Fig. 4(a)) has a higher reactants’ temperature and a lower flame temperature compared to OP2 (Fig. 4(b)). These conditions result in a higher temperature jump across the flame for OP2. The fuel mixture fraction distributions (not shown) at the burner mixing section exit (dump plane) are very similar for both OPs. In the OP2 case, more fuel is injected to reach the higher flame temperature. The normalized HRR on the centerplane is given in Fig. 4(e)for OP1 and in Fig. 4(f)for OP2. The HRR is normalized by the maximum value on the centerplane. OP1 shows a lower and more distributed HRR and a longer flame. Further, the radial temperature field at three different axial locations is shown in Fig. 4(c)for OP1 and Fig. 4(d) for OP2. The lower inlet temperature and higher exit temperature of OP2 are clearly visible. For completeness, the HRR cloud is given in Fig. 4(g)for OP1 and Fig. 4(h)for OP2. As is apparent from the centerplane HRR figures, the flame shapes are different. 3.2 Combustor Eigenmodes. The frequencies of thermoacoustic instabilities are usually close to the natural acoustic eigenfrequencies of the combustor [15]. The acoustic eigenmodes can be computed by solving the homogeneous Helmholtz equation (Eq. (4)) in a 3D FEM solver while considering a time-averaged temperature distribution of the flame within the combustor r� 1 q0r^p � �þx2 cp0 ^p¼0(4) The Helmholtz equation is obtained by Fourier transforming the homogeneous wave equation into the frequency domain (eixt Fig. 2 The probability density function of the pressure measurement PðgÞ, the time signal, and the corresponding frequency spectrum for (a) OP1 and (b) OP2 Fig. 3 (a) Computational domain (portion) of the CFD simulation and (b) longitudinal plane cut of the mesh showing its resolution together with the temperature contour plot Journal of Engineering for Gas Turbines and Power DECEMBER 2025, Vol. 147 / 121019-3 Downloaded from http://asmedigitalcollection.asme.org/gasturbinespower/article-pdf/147/12/121019/7530914/gtp-25-1285.pdf by guest on 06 October 2025 convention). The wave equation itself can be derived using the linearized massand momentum balance together with the equation of state. For the derivation, it is assumed that the flow velocity is negligible, i.e., low compared to the speed of sound. Further, any viscous losses as well as heat transfer are neglected. Ideal gas has been assumed, and gas properties (heat capacities cpand cV) are assumed constant. Small perturbations are postulated (linearity). Using the Helmholtz equation, the eigenmodes can be computed which are also referred to as the passive flame eigenmodes because the flame is only accounted for by incorporating the temperature field in the calculation and not as a HRR source term. The temperature distribution within the combustor is obtained from the CFD mean fields by area-weighted averaging across the combustor diameter, which results in a 1D temperature field along the x-axis for each OP. Here, we investigate five specific transverse eigenmodes which correspond to the first, second, and third azimuthal eigenmodes (1T, 2T, and 3T), as well as the first radial (1R) and first radial plus first longitudinal (1R1L) modes. They are shown in Fig. 5together with the eigenmode name abbreviations and the normalized frequencies. The acoustic boundary conditions correspond to a pressure release at the inlet ðp0¼0Þand fully choked ðu0¼0Þat the outlet. The investigated combustor eigenmodes (Fig. 5) exhibit a high degree of rotational symmetry which is shown in Fig. 6. Therefore, if the pressure along multiple cross sections (Fig. 6(a)) of a sampling plane (Fig. 6(c)) of the combustor (parallel to the symmetry axis) is normalized, the curves are close to identical (Fig. 6(b)), especially in the region close to the flame. Thus, the 3D pressure field can be expressed by a single normalized 1D modeshape which is locally scaled with a pressure value p0ðxref ,y,zÞat a reference plane perpendicular to the combustor axis. In this case, the reference plane is close to the dump plane. A single normalized pressure curve along a cross section can be computed for all five eigenmodes. This results in a total of five 1D pressure modeshapes, shown in Fig. 7within the combustor region between dump plane and outlet. These normalized eigenmode shapes can subsequently be used to compute the reheat flame’s HRR response to the acoustic pressure. 3.3 Lagrangian Framework for Auto-Ignition Flame Response Computation. To characterize the flame response to acoustic perturbations, flame transfer functions (FTFs) are used. FTFs relate the flame’s HRR response to an acoustic perturbation at a reference location upstream of the flame. The FTF is a frequency dependent complex valued quantity whose absolute value represents the gain in relation to the acoustic perturbation [16], i.e., how strongly the flame reacts to a given acoustic amplitude. The phase of the FTF describes the corresponding phase difference. Usually, the FTFs are given with respect to the acoustic velocity perturbation (u0) [17]. While this is applicable for planar wave perturbed propagationstabilized flames, which are sensitive to velocity perturbations, it was found that reheat flames are rather insensitive to velocity perturbations [18]. This can be seen by the low FTF gain for velocity perturbations shown by Gant et al. [19]. In comparison, the HRR gains in response to pressure and temperature perturbations are significantly higher [19]. The high sensitivity of auto-ignitionstabilized flames to temperature and pressure fluctuations is due to the sensitivity of the ignition delay time to small changes in temperature. Thus, due to the different sensitivity of auto-ignition flames, it is convenient to compute the FTFs with respect to the acoustic pressure p0instead of u0. Flame transfer functions can be deduced from experimental measurements where the flame is acoustically forced and the HRR response is measured [11]. Alternatively, CFD can be used [20]. Further, analytical and numerical models exist for certain flame configurations, for example, the FTF model from Zellhuber et al. [21], Gant et al. [19], and Gopalakrishnan et al. [13,22,23] that relate the HRR response of a 1D reheat flame to planar acoustic and entropic waves. To compute the FTFs for a transversely perturbed noncompact auto-ignition flame, the 1D-Lagrangian particle tracking framework from Gopalakrishnan et al. [13,23] is used in this work. This is done by simulating the ignition process using a set of 0D homogeneous constant pressure reactors in CANTERA. To compute the flame response for a given frequency, a set of particles are injected at the inlet at different time instances corresponding to the different phases Fig. 4 Temperature field on the centerplane for (a) OP1 and (b) OP2. Radial temperature field at three different axial locations: (c) OP1 and (d) OP2. HRR on the centerplane for (e) OP1 and (f) OP2. HRR cloud for (g) OP1 and (h) OP2. Fig. 5 Five passive flame eigenmodes for OP1 and their corresponding nomenclatures and frequencies 121019-4 / Vol. 147, DECEMBER 2025 Transactions of the ASME Downloaded from http://asmedigitalcollection.asme.org/gasturbinespower/article-pdf/147/12/121019/7530914/gtp-25-1285.pdf by guest on 06 October 2025 of the perturbing pressure oscillation. As the particles are transported downstream, they are perturbed by an acoustic wave with a locally varying pressure amplitude as defined by a 1D modeshape (Fig. 7). The governing equations of momentum, energy, and species mass balance for each injected particle in the presence of acoustic fluctuations are integrated in time to obtain the time evolution of the particle temperature. From this data, the ignition delay time as well as its modulation is deduced. The FTFs are computed by relating the pressure perturbation at the reference location (close to the dump plane) to the integrated HRR perturbations. The flame response is divided into two separate FTFs. FwXðxÞrepresents the gain and phase of the HRR (b_ qa,wX, Eq. (5)) for an eigenmode with an eigenfrequency Xand the corresponding modeshape w. GwXðxÞdescribes the frequency response of the axial movement (^xig,wX, Eq. (5)) of the flame. The FTFs are multiplied with the relative pressure perturbation ratio ^pðxrefÞ=p0to obtain the heat release rate and flame movement response. _ q0,ais the integrated HRR, and xig,0 is the mean ignition length. An analytical state space model of the computed FTFs is obtained via the VECTORFIT system identification tool [24]. This enables the flame response computation for complex frequency inputs and can subsequently be integrated in the finite-element solver b_ qa,wX _ q0,a¼FwXx ð Þ^p xref ð Þ p0 ^xig,wX xig,0 ¼GwXx ð Þ^p xref ð Þ p0 (5) The response of the auto-ignition flame to acoustic perturbations in a given geometrical configuration strongly depends on the nature of the acoustic field in the considered geometry. For a duct, the acoustic field can be written in terms of a simple mathematical expression. In more complicated geometries as considered in this work, the acoustic pressure field experienced by the particle cannot be written in terms of a simple mathematical expression. Due to this reason, the FTFs have to be computed on a mode-by-mode basis, where the acoustic pressure field introduced by each passive flame eigenmode is considered. This is then used in the Lagrangian framework to compute the particle evolution, and subsequently the heat release response of the flame is deduced. This is why the FTF expressions in Eq. (5) are valid for a specific eigenmode wX. The importance of modeling the flame response for each eigenmode separately was investigated in Ref. [14] for a 2D geometry. There it was shown, that if a pressure acts on the particle early on in the ignition process, it can heavily influence the ignition delay time. This consequently results in a large modulation of the ignition delay time and therefore in high FTF gains. In contrast, if the pressure only acts on the particle late in the ignition process, the ignition delay time only changes slightly, which leads to smaller FTF gains [14]. A comparison between two FTFs is shown in Fig. 8, where the influence of the 1T mode (diamonds) and a constant amplitude sinusoidal perturbation (circles) is compared. Differences in the flame response are visible in both the gain and the phase and can vary stronger for different eigenmodes or operating conditions. This underlines the importance of computing the FTFs for each mode separately under consideration of the acoustic eigenmode shape. To use the extended 1D-Lagrangian FTF framework, certain parameters are needed. The position of the flame in the 1Dframework is calculated from the 3D CFD mean field. Furthermore, the velocity is averaged between the dump plane and the mean flame position. During the 1D-Lagrangian framework computation, it is assumed that the pressure perturbation resulting from the eigenmode starts acting from the dump-plane onwards. This assumption is believed valid since high-frequency modes with transversal components are typically located locally in the combustion chamber, and the CFD progress variable of reaction remains very low in the mixing section. Figure 9shows the pressure perturbation (induced by the 1T modeshape) that the particle experiences during transportation with the mean flow. Differences are visible between the modeshape dependent amplitude (solid line) and a constant amplitude perturbation (dashed line), which results in the different FTFs of Fig. 8. The modeshape envelope is shown by the dashdotted line. Further assumptions are made in order to build the framework. First, the HRR perturbation due to acoustic velocity are neglected, and solely the influence of pressure, temperature, and density is accounted for (Eq. (6)) p0x,t ð Þ ¼AwXx ð Þ^ peixtþi/A T0x,t ð Þ ¼c−1 c T0 p0 p0x,t ð Þ q0x,t ð Þ ¼1 c2 0 p0x,t ð Þ (6) Fig. 6 Pressure of the 1T mode along different cross sections (a) for different points on the sampling cross section (c) on the dump plane. The pressure curves of (a) are normalized in (b). Fig. 7 Five normalized combustor eigenmodes for OP1 Fig. 8 FTF results from the Lagrangian framework for the 1T mode (diamonds) and a constant amplitude sinusoidal perturbation (circles) Journal of Engineering for Gas Turbines and Power DECEMBER 2025, Vol. 147 / 121019-5 Downloaded from http://asmedigitalcollection.asme.org/gasturbinespower/article-pdf/147/12/121019/7530914/gtp-25-1285.pdf by guest on 06 October 2025 The quantities p0ðx,tÞ, T0ðx,tÞ, and q0ðx,tÞcorrespond to the input values of the Lagrangian framework [13] and were adapted in Ref. [14] for the transverse case, which is also shown here. AwXðxÞ determines the location-dependent amplitude A of a certain modeshape wwith eigenfrequency X. The temperature and density perturbations result from this correspondingly [13]. Since velocity perturbations of purely transversal eigenmodes are perpendicular to the flow direction, it is assumed they have a negligible influence on the framework’s accuracy. Second, all modes are assumed to be standing and thus do not rotate or propagate. Given the boundary conditions of a pressure release inlet and choked outlet, as well as the investigated modeshapes, this is valid. Postprocessing of the pressure measurements showed that there is no rotating component for the investigated azimuthal modes. 3.4 Auto-Ignition Flame Response Modeling. In many cases, the flame can be assumed compact in longitudinal direction because the wavelength of planar acoustic waves (or longitudinal eigenmodes) is much larger than the flame thickness. Planar waves will acoustically perturb the flame uniformly, meaning that the acoustic perturbations are uniform across the flame front. In this case, the flame is acoustically compact. However, when the flame is transversely perturbed by the acoustic waves, the flame experiences a nonuniform pressure distribution. This means that different flame regions experience different acoustic perturbations, making the flame noncompact. To deal with this noncompactness, the flame is numerically divided into smaller subflames (blue dots in Fig. 10) in transverse zand y-directions. For simplicity, Fig. 10 only shows the segmentation in the z-direction, and not the entire yz-plane. This numerical segmentation allows the flame to be perturbed differently locally, depending on the modeshape. It is assumed that each subflame is perturbed by a specific value for the magnitude and phase of incident pressure wave (^ p). Consequently, for each subflame, a single value for the HRR perturbation (b_ qa,N) and axial location shift (^xig,N) with respect to the x-axis is computed. The segmentation is directly linked to the resolution of the numerical mesh, which, in this case, results in at least 40 elements per wavelength and multiple mesh elements in the flame brush. Each subflame is represented by a 1D Gaussian HRR distribution according to Eq. (7) and is defined for each yand z-axis combination with respect to the coordinate axis of the GT36 combustor (defined in Fig. 3(a)) _ q x,y,z ð Þ ¼_ q0,ay,z ð Þ r_ qy,z ð Þ ffiffiffiffiffiffi 2p pexp − 1 2 x−xig,0 y,z ð Þ r_ qy,z ð Þ !2 0 @1 A(7) The mean flame position xig,0ðy,zÞ, flame width r_ qðy,zÞ, and integrated HRR _ q0,aðy,zÞare calculated from the 3D CFD mean field for each operating point and are shown in Fig. 11. This is done by computing the parameters using a Gaussian fit (multiple studies used Gaussian distributions [14,19,21,23,25]) of the CFD volumetric HRR for each 1D cross section with constant yand z-coordinates in 3D, as was shown by Heinzmann et al. in the 2D case [14]. Thus, Fig. 11(a)indicates where on the yz-plane the stronger and weaker HRR regions are. Figure 11(b)defines the width of the HRR in x-direction for the 3D combustor volume. Figure 11(c)specifies the mean flame location xig,0ðy,zÞfor the 3D combustor volume. Using these three parameters, the 3D volumetric HRR zone in the combustor volume can be rebuilt analytically using Gaussian distributions in xdirection for every point on the yz-plane. Mathematically, the equation for a transversely noncompact autoignition flame that is perturbed by an acoustic combustor eigenmode is given by Eq. (12), whose formulation is explained in the following. The instantaneous HRR is obtained by accounting for the perturbations in the spatially integrated HRR b_ qa,wXy,z ð Þ ¼_ q0,ay,z ð ÞFwXx ð Þ^ p xref,y,z ð Þ p0 (8) along with the fluctuations in the flame motion ^xig,wXy,z ð Þ ¼xig,0 y,z ð ÞGwXx ð Þ^ p xref,y,z ð Þ p0 (9) and is spatially distributed using a Gaussian kernel as follows: Fig. 9 Pressure perturbation experienced by the particle while being transported by the mean flow. The perturbation induced by the 1T modeshape is plotted by the solid line and the constant amplitude perturbation in by the dashed line. The modeshape envelope is shown by the dash-dotted line. Fig. 10 Numerical segmentation method for a transversally perturbed flame in 2D 121019-6 / Vol. 147, DECEMBER 2025 Transactions of the ASME Downloaded from http://asmedigitalcollection.asme.org/gasturbinespower/article-pdf/147/12/121019/7530914/gtp-25-1285.pdf by guest on 06 October 2025 _ qwXx,y,z ð Þ ¼_ q0,ay,z ð Þ r_ qy,z ð Þ ffiffiffiffiffiffi 2p p… exp − 1 2 x−xig,0 y,z ð Þ−^xig,wXy,z ð Þ r_ qy,z ð Þ !2 0 @1 Ab_ qa,wXy,z ð Þ _ q0,ay,z ð Þ þ1 ! (10) The spatial distribution of the mean HRR is given by _ q0x,y,z ð Þ ¼_ q0,ay,z ð Þ r_ qy,z ð Þ ffiffiffiffiffiffi 2p pexp − 1 2 x−xig,0 y,z ð Þ r_ qy,z ð Þ !2 0 @1 A(11) The fluctuating HRR b_ qwXðx,y,zÞ, which is the source term for the Helmholtz equation, is consequently computed by subtracting the mean HRR from the instantaneous counterpart (Eq. (12)) b_ qwXðx,y,zÞ ¼ _ qwXðx,y,zÞ−_ q0ðx,y,zÞ(12) The final expression for the inhomogeneous Helmholtz equation for a 3D transversely perturbed and moving auto-ignition flame (Eq. (13)) is obtained by substituting the derived HRR source term (Eq. (12)) on the right-hand side of the Helmholtz equation (Eq. (4)). The equation can be numerically solved in a FE-solver in the frequency domain to extract the eigenfrequencies and linear growth rates. In this study, COMSOL MULTIPHYSICS 6.2 is used r� 1 q0r^ p � �þx2 cp0 ^ p¼−ixc−1 cb_ qwXx,y,z ð Þ (13) 4 Results In this section, we discuss the results obtained when using the framework developed in Sec. 3 to predict the linear stability of the GT36 CBB combustor. First, Sec. 4.1 presents the experimental measurements which serve as a benchmark to validate the model. Section 4.2 subsequently presents the obtained results from the framework. 4.1 Measurements of Pressure Pulsations in the GT36 Combustor. At both investigated OPs, distinct resonances can be observed. For the stable operating point, the resonances are clearly visible in the frequency spectrum of the measurement data, as shown in Fig. 12 measured by the two sensors 1 and 2. Sensor 1 is located at the center of the dump plane, while sensor 2 is located on the wall of the sequential liner close to the dump plane. The sensor locations are also indicated in Fig. 13(b). Figure 12 indicates that due to the fact that sensor 1 is located in the center of the dump plane, it is not possible to measure azimuthal modes. This is because the center of the dump plane coincides with the nodal lines in the pressure field of the 1T, 2T, and 3T modes (Fig. 5). In the contrary, sensor 2 on the sequential liner can measure the azimuthal modes’ pressure fluctuations, observable by the three peaks (dashed black line). If, however, a radial mode is present, sensor 1 can measure its pressure fluctuations, indicated by the two peaks (solid blue line) of the 1R and 1R1 L modes. Table 1compares the eigenfrequencies of the 1R and 1R1 L modes obtained from the FEM to the corresponding measurements of sensor 1. The predicted eigenfrequencies match very well with the measured ones. To ascertain the azimuthal and longitudinal nature of the predicted and measured eigenmodes and to ensure a comparison between modes of similar qualitative nature, it is instructive to analyze the phase relation between the pressure fluctuations at different locations. This is done by computing transfer functions between different sensors, both for the measurements and the FEM computation. Figure 13(a)shows the transfer function between sensor 1 in the very center of the dump plane and sensor 2 on the sequential liner wall at a small axial distance (indicated in Fig. 13(b)). The radius of the polar graph represents the ratio of the absolute pressure values, and the angle shows the phase difference between the fluctuations at the two considered locations. In the case of the FEM (*), a phase difference of 180 deg is present. This matches the theoretical expectation perfectly, as the antinodes of the 1R mode are located in the center and on the outer walls. A very close value of 200 deg is observed from the experimental data (�). The good match between the phase difference and in the eigenfrequency confirms that the measured resonance peak at 0.894 corresponds to the 1R mode. Subsequently, the transfer function is also computed for the mode with a scaled frequency of 1.033 to confirm the 1R1 L nature. Figure 13(c)shows the same analysis between sensor 1 and sensor 2. Identical phases are computed as for the frequency of the first radial mode, with the only difference in transfer function amplitude. This occurs because sensor 2 is located close to the nodal line in the FEM. For the investigated mode, a second transfer function (Fig. 13(d)) is computed between sensor 2 and sensor 3, which is also located on the sequential liner wall, but further downstream. Again, the FEM predicts the phase difference to be 180 deg, as is to be expected due to the antinode placement. Due to the close proximity of the sensor to the nodal line, the amplitude ratio is high. The transfer function phase of the experimental data however is 240 deg. This does not perfectly match the expected result of 180 deg, but clearly indicates that the pressures at these two locations are not in-phase. Hence, it is confirmed that the mode in question has a radial and longitudinal component present in its modeshape, and the excellent match in scaled frequencies confirms the 1R1 L nature. The 1R1 L mode is further analyzed using time domain analysis methods shown in Ref. [26]. Figure 14 displays the probability density function of the bandpass-filtered (60.037 around frequency peak) pressure measurement signal PðgÞ, the time signal, probability density function of the amplitude PðAÞ, the frequency spectrum, and the joint probability density function Pðg,_g=x0Þfor (a) OP1 and (b) OP2. All values are normalized with respect to the measured pressure at OP1. Clear differences are visible between the results of Fig. 12 Scaled frequency spectrum of the pressure measurements from sensor 1 located at the center of the dump plane (solid blue line) and from sensor 2 located on the wall of the sequential liner (dashed black line) for operating condition OP2line Fig. 11 Three flame characterization parameters: (a) integrated HRR _ q0,a(y,z), (b) flame width r_ q(y,z), and (c) mean flame location xig,0(y,z) for OP1 Journal of Engineering for Gas Turbines and Power DECEMBER 2025, Vol. 147 / 121019-7 Downloaded from http://asmedigitalcollection.asme.org/gasturbinespower/article-pdf/147/12/121019/7530914/gtp-25-1285.pdf by guest on 06 October 2025 OP1 and OP2. Whereas PðgÞis spread across a wider range of pressures for OP1, it is concentrated on a narrower region for OP2. Generally, for a thermoacoustic instability undergoing limit cycle oscillations, a bimodal shape would be expected for PðgÞ[27–30], which would show as a circular shaped ring in the g,_ g-plane [27,29,30]. This exact behavior cannot be observed. This confirms that the thermoacoustic instability is not at limit-cycle oscillation, as can also be seen by the intermittent nature in the time signal. Nonetheless, the high-amplitudes measured together with the broad PðgÞcurve and the non-Gaussian-like Pðg,_g=x0Þplot confirm the thermoacoustic instability. The time signal analysis for OP2 clearly indicates a stable mode, for which the acoustic resonances can be seen. The Gaussian-like shape in the ðg,_ gÞ-plane together with the narrow peak in the PðgÞconfirm this. 4.2 Linear Stability Analysis of Transverse Thermoacoustic Modes. The linear stability of the presented five eigenmodes is computed using the developed framework. The following points need to be kept in mind when interpreting the results: �The FEM is computed without any damping effects, neglecting viscous and thermal damping induced by the wall boundary layers. These are accounted for a posteriori by the method proposed by Romero et al. [31] when interpreting the FEM results. �In the GT36 combustor configuration considered here, specific dampers were installed at optimal locations to effectively dampen the 1T, 2T, and 3T eigenmodes. The dampers are not incorporated in the FEM computation, and their effect is considered a posteriori according to Ref. [32]. �An assessment of the damper performance [33] in terms of growth rate reduction shows that the 1R1 L mode is damped significantly less than the 1T, 2T, and 3T modes. �We compare linear growth rates from the experiment and the model below. We also show pulsation spectra for the two operating conditions. It should be noted that while the presented model is not able to predict limit-cycle amplitudes, the model can be qualitatively compared to pressure spectra to assess the stability of the system, even if experimental growth rates are not available. This allows to give reliable estimates of the thermoacoustic behavior helping in hardware modifications or defining the necessity of dampers with a very costeffective simulation. Fig. 13 Transfer function results between two pressure sensors at different locations, both for the measurement and the FEM computation. (a) Sensor 1 and sensor 2 for the 1R mode, (b) sensor configuration and legend, (c) sensor 1 and sensor 2 for the 1R1 L mode, and (d) sensor 2 and sensor 3 for the 1R1 L mode. Table 1 Comparison between the scaled resonance frequencies of the measurements and the FEM computation for OP2 Mode Scaled resonance frequencies Measurement FEM First radial 0.894 0.909 First radial þfirst longitudinal 1.033 1.049 Fig. 14 The pressure measurement time signal is postprocessed for one distinct frequency peak appearing in the unfiltered frequency spectrum of the time signal. The probability density function of the pressure measurements P(g), the time signal, probability density function of the amplitude P(A), the frequency spectrum, and the joint probability density function P(g,_ g=x0) are shown. All values are normalized with respect to the measurements at OP1. The plots show the analysis for the 1R1L mode for (a) OP1 and (b) OP2. 121019-8 / Vol. 147, DECEMBER 2025 Transactions of the ASME Downloaded from http://asmedigitalcollection.asme.org/gasturbinespower/article-pdf/147/12/121019/7530914/gtp-25-1285.pdf by guest on 06 October 2025 Figure 15 shows the stability prediction for (a) OP1 and (b) OP2. The small symbols on the zero growth rate line represent the passive flame eigenmodes. Analyzing the unstable operating condition OP1 first, it is visible that the 1T, 2T, and 3T eigenmodes are predicted as strongly unstable, which for the aforementioned reasons, is to be expected. This is in line with experience from previous measurement campaigns indicating that the 1T and 2T modes are typically more unstable than the investigated 1R1 L mode. Furthermore, the 1R mode is predicted to be stable and also does not show an instability in the measurements. Finally, the 1R1 L mode is predicted thermoacoustically unstable at a frequency of 1.005 and a linear growth rate of 1.31% (growth rate divided by its frequency of oscillation). The experimental growth rate is deduced by fitting an exponential curve to the bandpass-filtered pressure bursts [34], resulting in a mean value of 0.49%. If the growth rate reduction of the dampers is accounted for in the model as is described in Ref. [32] and the boundary layer losses are calculated by the method proposed in Ref. [31], the resulting net growth rate is reduced to 0.99%. The difference in growth rate of 0.50% is likely due to underestimated damping values and the simplifying assumption of our methodology. In essence, the framework is able to predict the thermoacoustic instability at a frequency of 1.005 as well as to identify its mode shape, confirmed by the experimental data in Sec. 4.1. According to the predicted growth rates of the stable OP2 in Fig. 15(b), the 1R1 L mode is stable. This is in line with the measurements (Fig. 14(b)), for which no instability is observed at that frequency for OP2. The 1R mode is predicted to be unstable. This mode is believed to be damped in the real combustor, as the dampers of the 3T mode are at ideal locations and at a similar damping frequency as the 1R mode. In summary, the framework is able to predict that the main pressure pulsation at a frequency of 1 is the 1R1 L mode. Furthermore, it shows that this mode is unstable at OP1 and stable at OP2. An additional analysis is performed, to investigate whether similar stability predictions considering only the shorter cylindrical section of the sequential liner just downstream of the area expansion can be made. This investigation is interesting as many of the investigated eigenmodes are particularly present in this cylindrical section and the smaller domain would result in a reduction of computational time. The computations for a shorter cylinder are done for OP1 and compared to the results of the real geometry. For this configuration, the inlet is characterized as a wall. The outlet is set to anechoic to match it as closely as possible to the full-length geometry. Using a pressure release or wall outlet boundary condition would place a pressure node or antinode at the outlet, which would not correspond to the condition of the full geometry. Figure 16 compares the stability predictions from (a) the symmetric GT36 can combustor and (b) the simplified geometry. While the 1T and 2T modes are captured well by the simpler geometry, the thermoacoustic behavior for the other modes cannot be reproduced. Thus, the downstream boundary condition of the full combustion chamber cannot be replaced by an artificial outlet further upstream with anechoic boundary condition. This is because the reflectiveness that the combustion chamber outlet imposes on the acoustic eigenmodes is not reproduced by the anechoic boundary condition. Therefore, the entire combustion chamber needs to be simulated, as was initially done in Fig. 15. 5 Conclusions In this paper, we present a 3D framework based on flame transfer functions and FEM to compute the eigenfrequencies and linear growth rates of transverse thermoacoustic eigenmodes in an industrial reheat combustion chamber at engine operating conditions. The thermoacoustic modeling approach captures the essential underlying physics of reheat flame dynamics and consists of multiple components: (1) a reactive CFD calculation performed in 3D to obtain information on the flame location, shape, the local HRR magnitudes, and the temperature distribution in the combustion chamber, (2) a Lagrangian particle-tracking algorithm to obtain the reheat flame response to acoustic perturbations, and (3) a numerical segmentation of the flame to capture the transversal noncompactness of the flame. The mathematical expression for a transversely noncompact reheat flame is then solved in a 3D Helmholtz solver. The framework is able to identify the unstable mode(s) and predict the linear growth rate, extracted from measurements, of an unstable mode with good accuracy. The validation is performed against experimental measurements for an unstable and a stable operating point. The presented method can be used to support the industrial development process of novel gas turbine reheat combustion chambers. In the early design phase, especially, it can be used to select from multiple hardware variants in a computationally very efficient way or to foresee dampers at relevant frequencies. Acknowledgment The research work presented in this paper was carried out within the FLEX4H2 project. The FLEX4H2 project was supported by the Clean Hydrogen Partnership and its members European Union, Hydrogen Europe, and Hydrogen Europe Research (GA Fig. 15 Stability predictions: (a) OP1 and (b) OP2. Small symbols represent the passive flame eigenmodes, and large symbols represent the active flame eigenmodes. The growth rate is given as a percentage of the respective active flame eigenfrequency. Fig. 16 Stability predictions: (a) GT36 can combustor and (b) short cylinder. Small symbols represent the passive flame eigenmodes, and large symbols represent the active flame eigenmodes. The growth rate is given as a percentage of the respective active flame eigenfrequency. Journal of Engineering for Gas Turbines and Power DECEMBER 2025, Vol. 147 / 121019-9 Downloaded from http://asmedigitalcollection.asme.org/gasturbinespower/article-pdf/147/12/121019/7530914/gtp-25-1285.pdf by guest on 06 October 2025