Abstract
Large eddy simulation data of a bluff-body stabilized flame are analyzed using spectral proper orthogonal decomposition (SPOD) to investigate: (i) the role of flame-vortex interactions in the dominant flow dynamics and (ii) how the proper choice of the cross-spectral density (CSD) defining SPOD can assist in identifying the underlying dynamics. Bluff-body flame holders aim to achieve stable flames under lean premixed conditions to minimize pollutant emissions. The recirculation region induced by the body promotes the mixing of hot combustion products with unburnt gases, preventing the global blowoff. However, the coupling between the shear layers and flame-induced vorticity sources can result in large flow structures that either contribute to increased flame stability or exhibit features typical of the early stages of flame blowout. SPOD is a data-driven technique remarkably powerful in extracting low-dimensional models. For each frequency, it computes a basis of orthogonal modes that maximizes the content of a predefined CSD in the leading modes. By choosing physically relevant variables to construct the CSD, different physics can be explored, which is used here to investigate the coupled dynamics between the flame-induced baroclinic torque, vortical structures, and the temperature field. The results show that the vorticity and temperature fields exhibit low-dimensional dynamics characterized by a narrowband frequency and its harmonics; these dynamics are varicose oscillations of the flame region, governed by the baroclinic torque. Sinuous oscillations typical of wake instability for nonreactive flows are also present, suggesting a competition between them.
1 Introduction
1.1 Dynamics of Reacting Wakes.
Bluff-body flame holders are of particular significance in the aerospace propulsion community, designed to achieve stable flames under lean premixed conditions, which is crucial for minimizing pollutant emissions. The bluff body induces a low-speed recirculation region in the presence of a higher central jet flow, promoting flame stability by sustaining the mixing of hot combustion products with unburnt gases. Flames anchored to bluff bodies have been extensively studied, both experimentally and numerically [1–3]. Numerous prior studies have demonstrated a tight, two-way coupling between the dynamics of vortical structures in the turbulent shear layers and the stable sustaining of the chemical reactions, thus influencing flame extinction and lean blowoff [4–7].
Under nonreacting conditions, the near wake region tends to develop sinuous motions due to the Bénard-von Kármán (BVK) instability. These motions give rise to the antisymmetric shedding of vortices and the near-periodic ejection of packets of recirculated fluid [8]. The observation of the same dynamics in reacting flows has been associated with initial stages of global blowout, as the release of hot reaction products from the near wake leads to failures in the re-ignition of the flame near its base [2]. For instance, Chaudhuri et al. [9] reported that intense antisymmetric vortex shedding during blowoff was responsible for the extinction of the flame within the shear layer, resulting in flame reaction being localized only within the near wake region. Furthermore, the organization of the shed vortices into coherent symmetric structures may even trigger thermoacoustic instabilities, either as a forced response [2,10] or as intrinsic instabilities [11].
In the majority of studies on flame-vortex interactions occurring in bluff-body stabilized flames, vortex shedding has been described as having a symmetric mode, an asymmetric mode, or presenting both, depending on the flow and operational conditions [10]. The reaction's heat release gives rise to temperature and density gradients across the channel that, in combination with the pressure gradient, produce a baroclinic torque. Owing to the qualitative and quantitative changes due to the baroclinic torque, this flame-induced vorticity source is ultimately capable of altering the dominant flame/wake oscillations from sinuous to varicose [12,13]. As opposed to the antisymmetric vortex shedding resulting from the sinuous BVK instability, the symmetric shedding of vortices in the latter case preserves the hot flame core, resulting in a more stable flame and reduced emissions.
Recently, Whitman et al. [14] investigated the effect of the streamwise pressure gradient in the near wake of three geometries of prism bluff bodies. In the same vein, Yalcinkaya and Gungor [15] used large eddy simulations (LES) to study the interaction of the longitudinal pressure gradient on the flame-vortex interactions of bluff-body stabilized flames under lean premixed conditions, close to lean blowoff. The U.S. Air Force Research Laboratory (AFRL) combustor geometry [3] was used as a reference, which considers a bluff body with an equilateral triangle shape placed at the symmetry plane of a parallel-walled channel. The base geometry was modified to alter the mean streamwise pressure gradient, introducing a convergent or divergent channel section downstream of the body. The results provided new insights into the complex interactions between the flame and the turbulent flow structures, highlighting the role of the baroclinic torque on flame stability, in agreement with previous studies in the literature.
Large eddy simulations like those by Yalcinkaya and Gungor [15] are today capable of reproducing with very high accuracy the complex interactions between hydrodynamic and chemical processes and predicting the combustion performance and operational envelope of a simple burner like the AFRL geometry. However, the raw LES data is often too complex to directly extract deep physical insights or develop low-dimensional models based on the dominant physical mechanisms. Data-driven analysis and modeling techniques then become powerful tools in extracting physical information that might be masked by the complexity of the LES. One such technique is the proper orthogonal decomposition (POD) [16–18] (also known as principal component analysis or Karhunen–Loeve decomposition). POD is based on constructing a matrix of covariance between the flow field fluctuations across the ensemble of empirical data (e.g., LES data) and calculating its eigenvalues and eigenfunctions. The set of eigenvalues and eigenfunctions (POD values and POD modes, respectively) constitute an orthogonal base for the exact projection of the original data, in which the statistical content of the first modes with respect to the covariance matrix is maximized. The first applications of POD in fluid mechanics aimed at identifying and extracting large-scale coherent structures existing in nonreacting incompressible turbulence. Gradually, the idea of POD has proven similarly powerful in other physical problems, like aeroacoustics, which combine turbulence with acoustic propagation (e.g., Refs. [19–21]).
In its simplest form, sometimes referred to as space-only POD, POD modes are “static” flow fields, representative of the data ensemble but regardless of temporal information. However, when the data corresponds to a stationary process and properly sampled time-resolved data are available, the properties of POD can be further exploited by performing a block-wise discrete Fourier transform from the time domain to the frequency domain, similarly to Welch's method for estimating periodograms [22]. Then, the temporal Fourier modes for the same frequency but different data blocks are used to build the equivalent covariance matrix (now referred to as the cross-spectral density matrix, CSD), in a manner analogous to standard POD. This frequency-domain variant is known today as spectral POD (SPOD) [23], but earlier applications are present in the literature, e.g., Ref. [24], going back to the original work by Lumley [17]. The main advantage of SPOD over standard POD lies in the incorporation of temporal information; the SPOD modes are not “static” flow fields but are associated with a well-defined frequency. The different SPOD modes for the same frequency are spatially orthogonal to each other, and temporally orthogonal to the SPOD modes for a different frequency. This property offers remarkably improved capabilities in identifying dominant physics and educating representative low-dimensional data-driven models.
As a note of caution, a related but different variant of POD exists in the literature, which is also named spectral POD [25], and consists of combining a space-only POD with a temporal filter. This homonymous technique has also been applied to the study of combustion instability (e.g., Ref. [26]), but differs from the (frequency-domain) SPOD used herein.
Space-only POD and SPOD allow for a flexibility in the definition of the covariance or CSD matrix that has not been fully exploited yet in many fields. Early works aiming at the extraction of coherent structures in incompressible, nonreactive turbulent flows employed a POD based on the product of velocity components, identical to the kinetic energy. Later, the analysis of aero-acoustic data provided examples of how the consideration of other metrics, closer to the root of the physical mechanisms, can result in remarkably more efficient low-dimensional models [27,28]. To the best of the authors' knowledge, only three works have employed SPOD in the analysis of reactive flow [29–31]. Shoji et al. [29] employed a CSD definition based on the planar laser induced chemiluminescence of OH; Brouzet et al. [30] used direct numerical simulation data to build a CSD approximating the linearized fluctuation energy, including the kinetic energy but also contributions from the pressure and entropy fluctuations; Casel et al. [31] considered a CSD based on the fluctuating kinetic energy alone, but including the density inhomogeneities of the mean flow.
The modal decomposition recovered by the different formulations of the proper orthogonal decomposition, and thus the insight that can be obtained from them, can be strongly dependent on the definition of the inner product. The shift from space-only POD to (frequency-domain) SPOD can be regarded as a change in the inner product definition. Similarly, Ek et al. [32] proposed a so-called permuted-POD in which the correlation involves the streamwise flow direction and time and compared its results to those of space-only POD. In an analogous manner, the use of different ad hoc inner products addressing flow variables of particular interest to the flame-vortex interactions, like the temperature fluctuations or the baroclinic torque, can provide different modal decomposition and insight on the physics at hand for the same flow field. This possibility, that has been successfully exploited in aero-acoustics (as mentioned in the paragraph above), has not been investigated for turbulent flows involving flame-vortex interactions; this is the objective of this work.
The rest of the paper is organized as follows. Section 2 briefly describes the LES of Yalcinkaya and Gungor [15] that are the departure point of this work. Section 3 explains the fundamentals of POD and SPOD methods and details the variants used. Section 4 presents results of the application of three different flow variables to construct the CSD matrix for the SPOD. The general characteristics of the SPOD results are described before addressing their physical interpretation for the problem at hand. Finally, Sec. 5 summarizes the findings of the SPOD analysis.
2 Large Eddy Simulation Database
This study departs from the large eddy simulations presented in Ref. [15]. The geometry is based on the AFRL combustor geometry [3], in which a bluff body is located on the symmetry axis of a plane channel. The bluff body geometry is an equilateral triangle with side length of mm, located 4D away from the inlet. The normalized dimensions of the combustor are 22D in channel length, 3D in height at the inlet, and 2D in width. Yalcinkaya and Gungor [15] studied modifications of the base geometry in which the channel walls diverge or converge with a constant angle past the bluff body base, in order to impose changes to the axial pressure gradient. Only one of these modified geometries is considered here, in which the walls have a divergence semi-angle of 1.5 deg.
Premixed air-propane with equivalence ratio of 0.65 enters through the inlet at 310 K temperature. The adiabatic flame temperature for this equivalence ratio is K. The inlet bulk velocity is m/s and the corresponding Mach number is 0.084. The Reynolds number based on and D is 39,000. Wave transmissive boundary conditions are imposed at the outlet and periodicity is imposed on the spanwise direction. No-slip adiabatic conditions are imposed at the walls of the bluff body. No-slip conditions are imposed at the top and bottom channel walls, with fixed temperature of 310 K and 600 K, respectively. A companion simulation considering nonreactive conditions is also performed, for comparison. In this simulation, air enters through the inflow and the temperature at both walls is set to 310 K. The rest of the boundary conditions are kept the same as for the reactive case.
The subgrid-scale stress tensor is modeled using the Smagorinsky model [33] and the subgrid-scale enthalpy and convective mass flux are calculated using an eddy viscosity approach. The partially stirred reactor combustion model is used for the chemistry–turbulence interaction. The skeletal propane-air reaction mechanism Z66 [34] that includes 24 species and 66 reactions is employed. The specific heat of each species is calculated using polynomial approximations based on the JANAF thermochemical tables. The perfect gas assumption is used to calculate density and the Lewis number is assumed as unity.
The numerical simulations employ the reactingFOAM solver of OpenFOAM-v7, using the PIMPLE algorithm with first-order Euler scheme and adaptive time-step for the temporal integration. Diffusion and convection terms are discretized using second-order schemes. Details on the formulation, methodology, and convergence studies can be found in the original Ref. [15].
Figure 1 shows the spanwise- and time-averaged streamwise velocity fields corresponding to the reactive-flow simulation and the companion simulation with nonreactive flow. Figure 2 shows instantaneous and spanwise- and time-averaged fields of the spanwise vorticity and temperature, respectively, for the reactive-flow simulation. Unless stated otherwise, dimensional quantities are used throughout this work. Note that only a reduced portion of the complete computational domain used in the LES is shown here, while the actual domain extends both upstream and downstream. The subdomain shown is the one used in the ensuing SPOD analyses.

Streamwise velocity (ms). Left: reactive flow. Right: nonreactive flow. The solid line corresponds to . Horizontal and vertical axes are and , respectively.

Reactive flow. Top: spanwise vorticity (s). Bottom: temperature T (K). Left: snapshot. Right: spanwise and time average.
3 Proper Orthogonal Decomposition, Spectral Proper Orthogonal Decomposition and Extended Spectral Proper Orthogonal Decomposition
3.1 Proper Orthogonal Decomposition.
Proper orthogonal decomposition is a data-driven methodology aimed at providing low-dimensional descriptions of the flow dynamics by identifying existing coherent flow structures in empirical data [16,18]. In this work, POD is applied to time-resolved data from large eddy simulations, that is written in compact form as , where the vector q comprises all the flowfield variables either involved in the LES calculation or computed a posteriori. The vector comprises the streamwise x, transverse y, and spanwise z coordinates, and t is time. In Figs. 1 and 2 and subsequent analogous contour plots, x is the horizontal direction, y is the vertical direction, and z is perpendicular to the image.
where superscript * denotes complex conjugate transpose. When the inner product is applied to scalar components of q, e.g., it defines their spatial covariance tensor: . The tensor is approximated using the (finite number) realizations of the random process , taking the form of a covariance matrix which, by definition, is symmetric and positive semidefinite. Accordingly, it has an eigenvalue decomposition in which:
The eigenvalues are all real and non-negative, and can be ranked as .
The eigenfunctions form an orthogonal set, and are customarily normalized so that , where the latter is the Kronecker delta.
The properties of the expansion and the eigenmodes of result in the following properties:
The eigenvalues sum to the total expected variance of the data under the inner product defined by tensor W: , where denotes expected value.
As the are ranked in descending order, a truncated expansion captures a larger fraction of the variance than any other orthogonal expansion of the same order.
The orthogonality of the eigenfunctions imply that the expansion coefficients are uncorrelated: .
These properties make POD a powerful tool in the eduction and analysis of patterns existing in the flowfield data, as it aims to extract flow structures (POD functions) of high spatial covariance over the data ensemble, i.e., spatially coherent and representative across the realizations. The choice of the inner product via the weight tensor allows to seek for coherent structures associated with flow variables of particular relevance. The most widely used definition considers the summed product of the three velocity components (the kinetic energy), thus seeking for POD modes that maximize their kinetic energy content. Alternative definitions can be more relevant depending on the physics at hand. In this work, three different covariance definitions are employed, aimed at the study of the flame-turbulence interactions.
3.2 Spectral Proper Orthogonal Decomposition.
Proper orthogonal decomposition functions become Fourier modes on the spatial directions in which the realizations are translationally invariant (homogeneous) [18]. Analogously, if the stochastic process is statistically stationary in time, the time-resolved data can be transformed to the frequency domain, as is customarily done in time-signal analysis when using estimated periodograms. The LES data analyzed in this work is both homogeneous in the spanwise direction z and stationary in time, which is exploited in the analysis.
where is the spanwise wavenumber and is the spanwise domain length. In practice, a direct Fourier transform is applied to the LES data, taking advantage of the evenly spaced discretization on the spanwise direction. The spanwise wavenumber is thus discretized as , , with and is the number of discretization points on the spanwise direction.
Following a procedure analogous to Welch's method for periodogram estimation [22], the discrete time series are divided into blocks, each of them with an equal number of snapshots . A Hanning window is used to reduce aliasing and a 50% overlap between the blocks is used to reduce the spectral leakage. The windowed discrete Fourier transform of each block b is then computed using a fast Fourier transform algorithm, to obtain . The discrete frequencies depend on the sampling time-step and through the frequency bin . The maximum discrete frequency is determined by Nyquist criterion, .
Spectral POD [17,23] exploits the stationary nature of the stochastic process to incorporate directly the temporal information contained in the Fourier modes; a POD problem is defined independently for each mode which is uncorrelated to the other Fourier modes.
with being the discrete Hanning window.
By construction, the SPOD modes at each are spatially orthonormal, , the identity matrix. Modes of different are not spatially orthonormal but they are orthogonal in the original space, as they remain Fourier modes.
so that they contain their estimated fluctuation amplitudes.
3.3 Choice of Inner Product and Extended Spectral Proper Orthogonal Decomposition.
The former derivation of POD and SPOD considered that the snapshots q consist of an arbitrary number of scalar fields (velocity, vorticity, temperature, and so on) and a generic weight tensor . The straightforward implementation would yield intractably large matrices and . In practice, the weight tensor is defined to involve the CSD of a reduced number of variables. In this work, we consider three different definitions of the weight tensor, leading to three different CSD matrices: the spanwise vorticity component , temperature , and the spanwise component of the baroclinic torque . In each case, only one scalar field (, or ) is involved in the calculation of . This reduces drastically the computational cost, as only the variable directly involved in the inner product is stored in the matrix .
For each definition of the CSD matrix, SPOD computes a different set of and , based on the coherent fluctuations of the scalar field used in the inner product. However, it is possible to compute, a posteriori, the other flowfield variables associated with each SPOD mode, that were not involved in the inner product. This is achieved, after solving the eigenvalue problem for the calculation of and , by extending matrix to incorporate all the variables of interest before applying Eq. (9). This is usually known as extended POD (or extended SPOD), and expands the utility of traditional POD in understanding complex flow phenomena [36]. For the three different CSD definitions used here, the scalar fields , and are monitored in this work. SPOD is also applied to the nonreactive simulation for the sake of comparison; in this case, only the CSD definition based on is applied.
4 Results and Discussion
To facilitate the SPOD analysis, the original LES data, which is stored as element-center values of a structured finite volume mesh, are interpolated to an overlapped Cartesian mesh. Points in the Cartesian mesh that lie outside of the LES domain (appearing in white in Fig. 1 and subsequent) are marked and disregarded in the analysis. Preliminary studies showed that a resolution of 201 points in the streamwise domain and 401 points in the wall-normal domain is sufficient to deliver converged SPOD spectra at the relevant frequency range.
An undersampling of the LES data on time was done to keep the volume of data tractable while retaining its statistical significance. The stored time-step s was chosen to obtain a Nyquist frequency Hz. A total of snapshots were stored, extending for a total dimensional time of 1 s, or 20 flow-through-times after an initial ten flow-through-times for wash-out. The number of samples per block was chosen to yield a frequency bin Hz. An overlap of 50% between blocks was applied, resulting in a total of blocks.
where is the density, P is the static pressure, and is the unit vector on the spanwise direction, is calculated at each snapshot in the finite volume mesh and then processed together with the other variables of interest. For the nonreactive flow data, only the spanwise velocity field is considered.
Spectral proper orthogonal decomposition analyses are performed for three different definitions of the inner product. Preliminary results, not shown here, show that the power spectral density of the fluctuations decay fast with increasing spanwise Fourier number for all frequencies. The results presented in what follows consider only (spanwise-averaged fluctuations) and the first spanwise mode . It is anticipated that strongly coherent flow dynamics exist, which are dominated by the modes.
4.1 Reactive Flow
4.1.1 Spectral Proper Orthogonal Decomposition Based on Spanwise Vorticity.
Figure 3 shows the spectra of SPOD values for the CSD based on the spanwise vorticity and for the spanwise wavenumbers and . SPOD values are shown by the grayscale lines, going from black (the leading mode) to almost white. The thicker top lines in the figure (red in the online version) show the total power spectral density of .
At each , the relative significance of the individual SPOD mode is quantified using two numbers. The first one is , which shows the relative contribution of mode i with respect to the following one. This number is always . If is a large number, the contribution of mode to the CSD is negligible compared to that of mode i. The second quantity is and is always contained between 0 and 1. Values of close to one indicate that the leading mode accounts for nearly all the CSD at the given ; the leading POD mode provides then a good description of the fluctuation structure and the dynamics are said to be of low rank. Table 1 summarizes these numbers for the three leading SPOD modes for all the cases analyzed in this work.
Comparison of SPOD eigenvalues for reactive and nonreactive flows
Mode | Mode | Mode | ||||||||
---|---|---|---|---|---|---|---|---|---|---|
Flow; CSD variable | Azimuthal mode | |||||||||
Reactive; , Hz | 1.34 × 102 | 0.83 | 8.33 | 1.60 × 101 | 0.1 | 2.93 | 5.48 | 0.03 | 4.96 | |
0.32 | 0.13 | 1.32 | 0.24 | 0.1 | 1.11 | 0.21 | 0.09 | 1.13 | ||
Reactive; T, Hz | 5.03 × 101 | 0.88 | 9.23 | 5.45 | 0.1 | 6.98 | 0.78 | 0.01 | 3.12 | |
0.09 | 0.33 | 1.3 | 0.07 | 0.25 | 2.73 | 0.03 | 0.09 | 1.21 | ||
Reactive; , Hz | 6.44 × 106 | 0.77 | 17.48 | 3.69 × 106 | 0.04 | 1.18 | 3.12 × 106 | 0.04 | 2.04 | |
7.71 × 105 | 0.08 | 1.05 | 7.34 × 105 | 0.08 | 1.08 | 6.81 × 105 | 0.07 | 1.04 | ||
Nonreactive; , Hz | 2.59 × 102 | 0.95 | 6.4 × 101 | 4 | 0.01 | 1.77 | 2.25 | 0.008 | 1.72 | |
1.02 | 0.17 | 1.55 | 0.66 | 0.1 | 1.24 | 0.53 | 0.09 | 1.13 |
Mode | Mode | Mode | ||||||||
---|---|---|---|---|---|---|---|---|---|---|
Flow; CSD variable | Azimuthal mode | |||||||||
Reactive; , Hz | 1.34 × 102 | 0.83 | 8.33 | 1.60 × 101 | 0.1 | 2.93 | 5.48 | 0.03 | 4.96 | |
0.32 | 0.13 | 1.32 | 0.24 | 0.1 | 1.11 | 0.21 | 0.09 | 1.13 | ||
Reactive; T, Hz | 5.03 × 101 | 0.88 | 9.23 | 5.45 | 0.1 | 6.98 | 0.78 | 0.01 | 3.12 | |
0.09 | 0.33 | 1.3 | 0.07 | 0.25 | 2.73 | 0.03 | 0.09 | 1.21 | ||
Reactive; , Hz | 6.44 × 106 | 0.77 | 17.48 | 3.69 × 106 | 0.04 | 1.18 | 3.12 × 106 | 0.04 | 2.04 | |
7.71 × 105 | 0.08 | 1.05 | 7.34 × 105 | 0.08 | 1.08 | 6.81 × 105 | 0.07 | 1.04 | ||
Nonreactive; , Hz | 2.59 × 102 | 0.95 | 6.4 × 101 | 4 | 0.01 | 1.77 | 2.25 | 0.008 | 1.72 | |
1.02 | 0.17 | 1.55 | 0.66 | 0.1 | 1.24 | 0.53 | 0.09 | 1.13 |
The spectra for shows dominant narrow-band frequencies; the fundamental frequency occurs at Hz ( Hz, owing to the frequency bar Hz) and a number of harmonics with decreasing amplitudes is observed. Such spectra is typical of flow dynamics dominated by a flow oscillator behavior, like the Bénard-von Kármán instability [37–39]: a self-sustained hydrodynamic instability gives rise to flow fluctuations at a well-defined frequency. Then, nonlinear interactions originate harmonics and a mean flow distortion, ultimately leading to the saturation of the fluctuation amplitude.
The spectra in Fig. 3 shows that the first SPOD mode accounts for nearly all the CSD at the narrow-band peaks. At these frequencies, is a large number and is close to 1. Figure 4 shows the leading SPOD mode for the fundamental frequency Hz. The real components of the scalar fields , and are shown. The modes resemble a wavy structure that propagates along the streamwise direction, presenting a sequence of maxima and minima. Note that the SPOD modes are complex quantities for nonzero frequencies. The imaginary components (not shown) are qualitatively similar to the real ones, but with their respective peaks being phase-shifted on x with respect to those of the real component. It should be kept in mind that the SPOD modes describe fluctuations of the variables over the mean field; negative values of T or are then physically meaningful. The fields and T for this mode are nearly exactly symmetric, while the corresponding fields and are approximately antisymmetric with respect to the bluff-body axis. This combination of symmetries corresponds to coupled varicose oscillations of the bluff-body wake and flame.
These symmetries are not imposed, but recovered implicity by the SPOD. Indeed, the vorticity and baroclinic torque fields present an evident lack of symmetric/antisymmetric behavior near the top and bottom channel walls. Other SPOD modes, shown below, present a similar lack of symmetry. These asymmetric structures partially stem from the lack of symmetry in the wall conditions imposed in the LES, in which the bottom wall is heated but the top wall is cold. The observation that SPOD modes are strongly symmetric or antisymmetric across except for the channel walls' boundary layers illustrates the ability of SPOD in extracting the coherent structures that dominate over the entire domain. Figure 5 compares the first three SPOD modes for the same Hz), which presents similar flow structures to the first SPOD mode but alternating the symmetric/antisymmetric structure.

Reactive flow. Leading SPOD modes based on spanwise vorticity for Hz and . Top: (s). Bottom: T (K).
The combination of symmetries and antisymmetries observed in the first POD mode (Fig. 4) constitutes a varicose oscillation of the wake/flame region: at a given streamwise coordinate x, the temperature fluctuation has the same sign on both sides of the axis. If this fluctuation is superimposed to the mean temperature field to produce a low-dimensional model of the time-dependent field (similar to Fig. 2), the isocontours of temperature would present a varicose structure, i.e., the cross-stream distance between the upper and lower flame fronts would increase and reduce successively along x. An analogous description applies to the other scalar flow variables related to the reaction, including , and the streamwise velocity component. Conversely, the corresponding transversal (on the y axis) velocity field would be antisymmetric, as would be the spanwise vorticity and its originating baroclinic torque. This velocity field corresponds to a succession of contractions and expansions of the wake along the streamwise direction. On the other hand, the combination of symmetries and antisymmetries observed in the second POD mode (Fig. 5, center panels) constitutes a sinuous oscillation of the wake/flame region; the bluff-body wake and the location of the flame fronts displace successively toward positive and negative y locations. This is consistent with antisymmetric fluctuations of the state variables and other scalars associated with species and chemical reactions, and symmetric spanwise vorticity fields.
Figure 6 is analogous to Fig. 5 but for the first harmonic frequency Hz and . The dynamics for these modes are again low-dimensional. The leading (and subsequent odd number) SPOD mode structure describes varicose motions of the flame, characterized by symmetry and T and antisymmetric and fields, while the even-number SPOD modes correspond to sinuous flame oscillations. The apparent wave-like structure of these modes presents a wavelength which is half of that for Hz. This is a typical result in the description in the frequency domain of coherent fluctuations that are convected by the flow: convective waves present a streamwise phase velocity which is weakly dependent on the frequency, , where is the axial wavenumber. Thus, increasingly high frequencies correspond to smaller wavelengths. The leading SPOD modes for the higher harmonics (not shown) recover this trend.

Reactive flow. Leading SPOD modes based on spanwise vorticity for Hz and . Top: (s). Bottom: T (K).
Note that the SPOD modes shown in the figures are scaled with their spectral amplitude, as described by Eq. (9). This, together with the color legends given in the figures, provides a visual comparison of the relative weight of each mode on the complete fluctuating field (e.g., Figs. 5 and 6). Observation of the SPOD modes reveals that, for the narrow-band peaks, the odd SPOD modes () correspond to varicose flame oscillations and the even modes () to sinuous oscillations. However, it is remarked that the leading mode, describing varicose oscillations, dominates the dynamics associated with the inner product. As will be shown later, the same varicose mode is covered for different inner products.
As opposed to the modes, the higher spanwise Fourier modes do not present low-order dynamics. As shown by Fig. 3 (right), the leading SPOD mode does not account for a substantial part of the CSD for any frequency. This also implies that the modes are not expected to present a well-defined spatial structure. Figure 7 shows the first SPOD at Hz for the spanwise mode .
4.1.2 Spectral Proper Orthogonal Decomposition Based on Temperature.
The previous subsection presented SPOD results based on the CSD of spanwise vorticity. The extended SPOD modes contain all the flowfield variables but the decomposition attends to coherence in the spanwise vorticity alone. While a strong coupling between the flame (chemical reactions, heat release, and temperature increase) and the turbulence (strongly nonlinear vortical motions involving a wide range of temporal and spatial scales) is expected, SPOD based on a different CSD can potentially recover a different low-dimensional description of the flow. This subsection presents results for a CSD based on the temperature fluctuations.
Figure 8 shows the SPOD values for the and spanwise modes. These spectra are qualitatively very similar to the SPOD spectra based on vorticity (cf. Fig. 3) but present quantitative differences. The absolute values of the and the PSD are not to be compared for different inner products, as they involve different magnitudes of the flow field. Two relevant aspects are observed: (i) while the dominant narrow-band peaks are again present for , the leading SPOD mode accounts for a sensibly smaller part of the CSD. This denotes a higher relevance (on the temperature field) of the second SPOD modes, which describe sinuous flame oscillations. (ii) As opposed to the results for the -based CSD, the combination of the SPOD modes 1 and 2 for modes now accounts for most of the CSD. This implies that the spanwise fluctuations of the temperature field are more coherent than those of the vorticity field, but they remain of considerably lower amplitude than the spanwise homogeneous ones.
The SPOD modes for the same pairs shown in the previous subsection are shown for the temperature-based CSD in Figs. 9–11, for comparison. For modes presenting a large , a strong resemblance is observed for all the variables (cf. Figs. 4 and 9). This suggests that the fluctuations of all the flowfield variables are strongly coupled and thus the coherence of, e.g., temperature fluctuations is accompanied by vorticity fluctuations and vice-versa. However, this does not necessarily imply a causal relation between them.
4.1.3 Spectral Proper Orthogonal Decomposition Based on Baroclinic Torque.
The previous two sections considered SPOD results based on the CSD of two variables ( and T) that are representative in the description of turbulence and flame dynamics. Owing to the coupled nature of flame-turbulence interactions, both analyses recover a very similar picture: varicose motions of the near-wake region, coincident with the flame front, with characteristic dimensions comparable to the bluff-body dimension. However, in the equivalent nonreacting flow, BVK instability dominates the vorticity dynamics, leading to sinuous (rather than varicose) motions, characterized symmetric vorticity field (as discussed in Sec. 4.2). The combustion reaction introduces two vorticity production mechanisms, namely, the baroclinic torque and the vortex dilatation. Of special interest here is the baroclinic torque that is generated along the entire flame front with different signs due to the density gradient between unburned reactants and products, and ultimately favors an antisymmetric vorticity field. In consequence, it is expected that the dominant varicose oscillations identified in the previous analyses stem from fluctuations of the baroclinic torque, and hence very low-order dynamics (equivalently, very high coherence) exist for the latter.
Figure 12 shows the SPOD values for the and spanwise modes, together with the PSD of the spanwise component of the baroclinic torque . These spectra present some relevant qualitative differences with respect their counterparts based on or T: (i) Except for the narrow-band peaks, the PSD does not decay appreciably across the frequency range analyzed. (ii) Similarly, except for the narrow-band peaks, SPOD does not recover a low-rank behavior as many SPOD modes are required to add to a quantity of the same order of the PSD. This is remarkable for the modes. (iii) At the narrow-band peak frequencies, the first SPOD mode accounts for nearly all the CSD, rendering the higher SPOD modes irrelevant. As expected, the leading SPOD mode at the peak frequencies depicts varicose flame motions (see Fig. 13). Because the flowfields of , T, and present such a low-rank behavior at these frequencies, their first SPOD modes are visually indistinguishable.

Reactive flow. SPOD spectra for CSD matrix based on the spanwise component of the baroclinic torque. Left: . Right: .

Reactive flow. First SPOD mode based on the spanwise component of the baroclinic torque at Hz and , , .
4.2 Nonreactive Flow: Spectral Proper Orthogonal Decomposition Based on Spanwise Vorticity.
For the sake of comparison, this section presents analogous SPOD results for LES data corresponding to nonreactive flow under the same incoming velocity conditions. Figure 14 shows the SPOD spectra for the and spanwise modes. Comparison of these spectra with those for the reactive-flow simulation (Fig. 3) leads to the following observations: (i) For , the nonreactive flow presents a faster decay of the spanwise vorticity PSD with frequency. A series of narrow-band peaks also appears; the fundamental frequency is reduced to Hz. This frequency is in agreement with previous experimental and numerical works on the Volvo burner geometry [40,41]. The amplitudes of the narrow-band peaks also decay with frequency faster than for the reactive case. (ii) For , the first SPOD modes at the narrow-band peak frequencies account for nearly all the PSD; e.g., as shown in Table 1, for the leading SPOD mode based on spanwise vorticity at the fundamental frequency for the reactive flow, while for the nonreactive flow increases up to . This renders the second and later SPOD modes as negligible with respect to the first one. (iii) The spectra for does not have evidence of narrow-band peaks, neither in the PSD or in the individual SPOD modes. The proximity of the SPOD values indicates the absence of spanwise inhomogeneous large scale coherent structures.
Figure 15 shows the first three SPOD modes for the fundamental frequency Hz and its first harmonic Hz for the nonreactive flow, to be compared to the corresponding reactive-flow counterparts in Figs. 5 and 6. Here again, some differences are observed. For the reactive flow, odd-number SPOD modes correspond to varicose flame oscillations and thus to spanwise vorticity fields that are antisymmetric about , while the even-number SPOD modes correspond to symmetric vorticity fluctuations. This was satisfied by the SPOD modes for both the fundamental frequency and its harmonic. Conversely, for the nonreactive flow, all the leading SPOD modes for the fundamental frequency present a symmetric vorticity field, while those for the harmonic frequencies present an antisymmetric one. Recalling that only the leading SPOD is relevant in the nonreactive case, the SPOD results show that the fundamental frequency consists of sinuous oscillations of the wake, typical of the BVK instability, while its first harmonic corresponds to a varicose oscillation.

Nonreactive flow: leading SPOD modes based on spanwise vorticity for Hz (top) and 240 Hz (bottom), and . The contours correspond to (s).
5 Conclusion
This paper studied the interaction between the flame dynamics and vorticity in a premixed flame anchored to a bluff-body. LES data were analyzed using spectral proper orthogonal decomposition, which provides a powerful methodology to identify and extract coherent structures. “Coherent” stands here for a high cross-spectral density of a predefined spatial inner product. Ad hoc inner products are tested in the analysis of the physics, aimed at identifying low-rank dynamics associated with different physical phenomena: spanwise vorticity as representative of turbulence, temperature as representative of the combustion, and the baroclinic torque as the underlying mechanism responsible for the selection of the preferred coupled oscillations. For the sake of comparison, SPOD is also applied to a companion simulation considering nonreactive flow at the same conditions.
The results for the reactive-flow case show a series of intense narrow-band peaks characterized by a fundamental frequency Hz and its harmonics. In terms of the Strouhal number defined with the base length of the bluff body D and the inflow velocity , this corresponds to . The fluctuation field is remarkably low-rank at these frequencies, with the leading SPOD mode accounting for over of the CSD () and being well separated from the second one ( for all peaks and inner products, reaching for the CSD based on the baroclinic torque). The leading SPOD mode describes varicose oscillations of the flame, which correspond to symmetric fluctuations of the temperature and OH mass fraction fields—representative of the chemical reaction—and antisymmetric fluctuations of the spanwise vorticity and baroclinic torque. The same qualitative description holds for the leading SPOD mode for the first harmonic frequency, Hz (). However, the second SPOD modes for the narrow-band peaks describe sinuous flame oscillations. The co-existence and competence between sinuous and varicose oscillation modes for reactive wake flows has been documented in the past (e.g., Ref. [10]), as is the dominance of the varicose oscillations recovered herein.
The results for the nonreactive flow show the same series of narrow-band peaks, but the fundamental frequency is reduced to Hz (). For the narrow-band peaks, the first SPOD mode virtually captures all the fluctuation PSD. Incidentally, and opposed to the reactive-flow results in which both sinuous and varicose oscillations were recovered by the successive SPOD modes for each frequency, the SPOD modes for the fundamental frequency recover sinuous wake oscillations and the SPOD modes for the first harmonic recover varicose wake oscillations. This suggests that, for nonreactive bluff-body flows, the wake oscillations are initiated by a sinuous global instability at the fundamental frequency, and its subsequent nonlinearities promote varicose oscillations at the first harmonic. Conversely, the coherent dynamics for the reactive flow are more complex and involve both sinuous and varicose oscillations at all narrow-band peaks. Present results do not provide insight on the nonlinear interactions between the different symmetries and frequencies.
The study provides evidence on the crucial role of the baroclinic torque in determining the preferred fluctuations that dominate the large scales of the turbulent flow and explain the global coupling between the flame and turbulent dynamics. Under nonreacting conditions, the Bénard-von Kármán instability dominates the turbulent flow dynamics, giving rise to sinuous fluctuations of the near-wake. Conversely, under the present reactive conditions, the heat released by the combustion establishes a density gradient which is symmetric across the bluff body. In combination with the streamwise pressure gradient along the channel, the density gradient gives rise to an overall antisymmetric baroclinic torque that induces large-scale antisymmetric spanwise vorticity, i.e., varicose motions of the near wake. This kind of oscillations is preferred toward flame stabilization on bluff-body anchored flames, as they favor the recirculation of hot reaction products toward the flame anchoring point, preventing global blowoff. However, while the varicose wake/flame oscillations are found here to be dominant, present results show that both kinds of oscillations co-exist. Furthermore, different studies in the literature suggest that even small variations of the reaction temperature [13], combustion chamber geometry [14,15], or the presence of additional acoustic disturbances [10,42] can promote either kind of oscillation.
Finally, a note of caution is given regarding the interpretation of the SPOD results. The low-rank description provided by SPOD does not imply that the whole of the flowfield fluctuations is coherent and “smooth.” The present flow is fully turbulent, and as such, three-dimensional and chaotic. The low-rank description should be interpreted as a model of the mechanism(s) from which the large-scale organized motions stem; additional dynamic mechanisms including vortex stretching, tilting and breakdown or thermodiffusive instabilities are also present and relevant to the flame dynamics, but their impact on the selection of the dominant flow structures is secondary. In this sense, our results illustrate how the choice of an inner product that describes more closely the physical mechanisms driving the coherent flow structures—the fluctuations of the baroclinic torque or temperature—is found to result in a significant reduction in the dimension of the low-rank description.
Funding Data
MIXSHY (Ref: TED2021-129719B-C21) funded by MCIN/AEI/10.13039/501100011033 and European Union's NextGenerationEU/PRTR (Funder ID: 10.13039/501100011033).
SHYGAS (Ref: PID2021-125812OB-C22) funded by MCIN/AEI/10.13039/501100011033 and European Union's FEDER (Funder ID: 10.13039/501100011033).
Horizon Europe MSCA postdoctoral fellowship (GA 101063992; Funder ID: 10.13039/501100000780).
Scientific and Technological Research Council of Turkey (TUBITAK) (Grant No. 219M139; Funder ID: 10.13039/501100004410).
Data Availability Statement
The datasets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request.
Nomenclature
- =
expansion coefficient of realization i in POD
- =
spanwise component of the baroclinic torque
- =
SPOD correlation matrix
- f =
frequency (Hz)
- =
correction factor in PSD calculations
- =
spanwise domain length
- =
number of blocks
- =
number of samples per block
- =
number of snapshots
- =
number of discretization points on the spanwise direction
- P =
static pressure (Pa)
- =
Fourier transform of flow variables matrix Q
- T =
temperature (K)
- =
spatial weight function
- =
OH mass fraction