Investigations and Improvement of Robustness of Reduced-Order Models of Reacting Flow
IInvestigations and Improvement of Robustness ofReduced-Order Models of Reacting Flow
Cheng Huang and Karthik Duraisamy University of Michigan, Ann Arbor, MI, 48109
Charles L. Merkle Purdue University, West Lafayette, IN, 47907
The impact of chemical reactions on the robustness and accuracy of projection-based Reduced-Order Models (ROMs) of fluid flows is investigated. Both Galerkinand Least-Squares Petrov Galerkin ROMs are shown to be less robust in reactingflows as compared to non-reacting flows. In particular, reacting flow ROMs show astrong sensitivity to the resolution and are often unstable. To identify the main under-lying causes, a representative problem that contains the essential physics encounteredin typical combustion dynamics problems is chosen. Comparisons with non-reactingsolutions are used to assess the impact of reactions. Investigations are focused onthree potential areas of significance: 1) preservation of conservation laws; 2) loss ofdissipation; and 3) existence of unphysical local phenomena. Results indicate thatconservation is relatively well-controlled and the global dissipation in the ROMs is ac-tually larger than that in the underlying CFD solutions. Spurious local phenomena are,however, highly deleterious. Specifically, the steep temperature gradients that char-acterize combustion can lead to oscillations in local temperatures even in the absenceof reactions. Representative calculations with physics-based temperature constraintsverify that eliminating such excursions results in considerable improvement in bothstability and future-state prediction capability. Assistant Research Scientist, Aerospace Engineering, [email protected]. Associate Professor, Aerospace Engineering, [email protected]. Professor Emeritus, School of Aeronautics and Astronautics, [email protected]. a r X i v : . [ phy s i c s . c o m p - ph ] J u l omenclature FOM = Full Order ModelROM = Reduced Order ModelPOD = Proper Orthogonal DecompositionLSPG = Least Squares Petrov Galerkin Q p = spatially discrete solution variables ˆR = spatially discrete residuals of full order model equations Ω ref ( t ) = boundary conditions ¯( · ) = time-averaged solutions ( · ) (cid:48) = mean subtracted unsteady solutions σ n = the n th singular value a n ( t ) = the n th POD temporal mode Φ n = the n th POD eigen-basis (cid:104)· , ·(cid:105) P = inner product operator f = spatially discrete residuals of reduced order model equations W n = the n th test function for Least Squares Petrov-Galerkin method A = pre-exponential factor k = kinetic energy s = entropyAD = artificial dissipation I. Introduction
Computational Fluid Dynamic (CFD) modeling of reacting flow serves an important role inthe investigation and understanding of combustion dynamics in aerospace propulsion systems op-erating at high pressure and temperature, the details of which cannot be quantitatively accessedthrough experiment. Accurate modeling of combustion dynamics in practical combustor designis key to improving engine performance, reducing failures and avoiding the devastating phenom-2na of combustion instability. Though modern computational capability has moved beyond theempirically-based design analyses of the past, high-fidelity simulations of full-scale practical com-bustors remain out of reach for day-to-day engineering design applications. While complete LESsimulations of combustion dynamics in laboratory-scale combustors are routinely pursued, [1], andrepresentative computations in small-scale full combustors have been reported [2–4], they remainprohibitively expensive. Therefore, for simulation-based design to be a reality, accurate and efficientmodeling methodologies are imperative.In the present work, we examine the use of model reduction techniques to achieve efficientand accurate modeling of combustion dynamics in complex geometries. Specifically, we pursueprojection-based Reduced Order Models (ROMs) [5–7] that couple Proper Orthogonal Decomposi-tion [8] (POD) with Galerkin and Petrov-Galerkin projection methods. Reduced-order models haveproven to be efficient in reducing complex partial differential equations (PDE) to low dimensionalordinary differential equations (ODE) in non-reacting flow problems such as flow control [9–11] andaeroelasticity [12, 13]. Recent studies have extended ROMs to reacting flow problems [14, 15] andpreliminary explorations of the POD/Galerkin technique have been carried out [16, 17] using a sim-ple 1D solver to establish a basic approach for ROM construction and its characteristics in reactingflows. More recently, the current authors initiated the first attempt to apply ROM techniques onrelevant combustor flow simulations [18, 19], which identified steep temperature gradients as one ofthe underlying challenges of reacting flow ROM development.This work is part of a larger initiative in which the authors are developing a Multi-FidelityFramework using Reduced Order Modeling (ROM) to enable efficient modeling of combustion dy-namics in practical full-scale engines. Using a rocket engine as an example, it should be noted thatcomplex flow phenomena are limited to the region in and around the injector element, where thepropellants mix and react through actions of shear, swirl, and/or direct impingement. In large,high-performance rocket engine combustors, where the problem of instability is particularly severe,hundreds of injector elements are typically used. The design and operation of these individual ele-ments control the spatio-temporal aspects of mixing, reaction, and the overall heat addition field.Although the number of injector elements is large (which precludes application of high-fidelity meth-3ds to the entire combustor), the elements themselves are rather small (the largest are on the orderof 1 cm). This fact can be used to an advantage by restricting the detailed characterization of theunsteady flowfield and its response to flow perturbations to one or a small set of injector elements.Then, a generalized reduced order model (ROM) for injector elements can be extracted from thesimulations and implemented in multi-fidelity design framework suitable for engineering analysis.The idea of such a multi-Fidelity framework has been successfully demonstrated in simpler modelproblems [20, 21].It is well-recognized that one of the most critical issues in ROM development is the lack of ro-bustness and stability when reducing a complex multi-scale system to a low-order ODE. As reported,the issues can come from the inherent lack of numerical stability in the POD/Galerkin method it-self [22], truncation of low-energy dissipative POD modes [23] and simplifications of higher-orderequations [24]. The balanced POD technique has been proposed to build numerically stable ROMfor linear systems [25, 26]. Bergmann et al. [23] proposed to add residuals of the Navier-Stokesequations to account for the absence of low-energy dissipative POD modes. Moreover, Lucia etal. [7] demonstrated the effectiveness of constructing stable ROMs by including additional artificialdissipation terms. Researchers have also tried to resolve the issues arising from the numerical prop-erties of the system equations. Rowley et al. [27] pointed out that defining a proper inner productcan be important when dealing with model reduction of the Navier-Stokes equations. Barone etal. [28–30] proposed to stabilize the reduced system by symmetrizing the higher-order PDE witha preconditioning matrix. For aeroelastic applications, Amsallen and Farhat [31] have shown theadvantages of using the descriptor form over the non-descriptor form of the governing equations.Carlberg et al. [32] have also demonstrated that the Least Squares Petrov-Galerkin (LSPG) methodgenerates ROMs that are consistently more stable than those generated by the Galerkin method viasymmetrization of the discrete Jacobian.All of the above research investigations were focused on the robustness and stability of ROMsfor non-reacting flow applications. The additional physical complexity arising from chemical reac-tions, however, further exacerbates the robustness issues that are observed in reduced-order models.Therefore, based on the challenges of previous reacting flow ROM development [18, 19], the main4bjectives of the present paper are to:1. Compare ROMs for reacting and non-reacting flow conditions based upon a test problem thatis representative of combustion dynamics problems;2. Reveal and narrow down the main underlying causes of ROM stability issues in reacting flowsimulations; and3. Implement and evaluate potential resolutions to improve ROM characteristics.Although our focus is on reacting flow problems, it is anticipated that observations gained from thepresent work will be useful for improving ROMs for non-reacting flows also.The remainder of the paper is organized as follows. In Section II, we briefly review the CFDmodel that is used to generate the datasets for ROM development and introduce the procedure forPOD mode generation and the Galerkin and Least Squares Petrov Galerkin (LSPG) model reductiontechniques. In Section III, we present the computational setup of the test problem used for theROM evaluations. In Section IV, we present computational results and assess the accuracy androbustness of ROMs developed from the benchmark test problem followed by a detailed diagnosis ofthe impact of combustion on ROM stability along with the implementation of a simple modificationthat provides immediate relief from the most severe symptoms. In the last section, we provideconcluding remarks and address future work.
II. FormulationA. Governing Equations
The computational infrastructure used for the full and reduced-order models solves conservationequations for mass, momentum, energy and species mass fractions in a coupled fashion and is basedon an in-house CFD code that has previously been used to model combustion instabilities [1]. Adetailed description of the Full Order Model (FOM) equations can be found in Appendix A. Thesemi-discretized version of the governing equations in terms of primitive variables yields the ODEsystem, ∂ Q p ∂t = ˆ R ( Q p , Ω ref ( t )) , (1)5here Q p ( t ) = [ Q p, ( t ) . . .Q p,i ( t ) . . .Q p,NI ( t )] T with N I as the total number of grid points and Q p,i ( t ) = [ p i ( t ) , u i ( t ) , T i ( t ) , Y k ( t )] T ∀ i . Equation 1 is a high dimensional ODE composed of( N I × N var ) equations where N var is the number of solution variables in Q p,i . Boundary condi-tions are applied at boundary faces in the finite volume scheme by the explicitly specified inputvector, Ω ref ( t ) . Moreover, it should be pointed out that both the FOM and ROM employ consis-tent finite volume discretization schemes and specifically for the current studies, second-order Roescheme [33] with flux limiter due to Barth [34] is used. B. Construction of POD bases for vector equations
POD bases are derived from a database created by storing snapshots of Full-Order Model (FOM)solutions of Eq. 1 over the entire computational domain at specific time instants. To enable theROM to faithfully predict the combustion response to various disturbances, multiple databases arecreated from FOM solutions, each of which is forced by driving a single (upstream or downstream)boundary condition, Ω ref ( t ) in a periodic manner while the remaining boundary conditions are heldfixed. Nonlinear effects are included by forcing at finite amplitudes.Given this database from FOM solutions of Eq 1, POD bases are calculated based upon thevector-valued method using the unsteady part of the variables, Q (cid:48) p ( t ) = Q p ( t ) − Q p (time-averagedvalues extracted before POD eigen-bases calculation), Q (cid:48) p ( t ) ≈ N p (cid:88) n =1 ˆ a n ( t ) σ n Φ n = N p (cid:88) n =1 a n ( t ) Φ n . (2)Here, σ n is the singular value and a n ( t ) the temporal amplitude of the n th POD mode of theorthonormal vector function, Φ n , where, (cid:104) Φ k , Φ n (cid:105) P = , if k = n , otherwise (3)with an inner product defined as (cid:104) u , v (cid:105) P = u T P v . A normalization matrix, P , must be applied toscale the vector-valued POD bases. Here, we normalize all fluctuation quantities by their maximumamplitude, P = diag (cid:16) ˆ P . . . ˆ P i . . . ˆ P NI (cid:17) , (4)6here ˆ P i = diag (cid:16) p (cid:48) max , u (cid:48) j,max , T (cid:48) max , Y (cid:48) l,max (cid:17) and ϕ (cid:48) max ≡ Max { ϕ (cid:48) ( x, t ) } , ∀ x min ≤ x ≤ x max and t > . C. Model Reduction
In this work, projection-based model reduction is accomplished using both Galerkin- and Least-Squares Petrov Galerkin-POD methods. Galerkin-based model reduction procedures are formulatedupon the continuous-time representation of Eq. 1, which is then projected onto the k th test function, V k , chosen to be the k th basis, Φ k , obtained in Eq. 2, (cid:28) V k , ∂ Q p ∂t (cid:29) P = (cid:68) V k , ˆ R ( Q p , Ω ref ( t )) (cid:69) P . (5)where the inner product operator is introduced in Eq. 3, and Q p is approximated using Eq. 2,following which a time-variant ODE system is obtained, d a ( t ) dt = f ( a ( t ) , Ω ref ( t )) , (6)where a ( t ) = (cid:2) a ( t ) · · · a i ( t ) · · · a N p ( t ) (cid:3) T , f ( t ) = (cid:2) f ( t ) · · · f i ( t ) · · · f N p ( t ) (cid:3) T , the dimension of theROM ODE in Eq. 6 is N p , orders of magnitudes smaller than the N I × N var degrees of freedomin Eq. 1. As indicated, the same boundary condition, Ω ref ( t ) , appears in the ROM ODE as in theCFD equations, Eq. 1.The Least Squares Petrov-Galerkin (LSPG) model-reduction technique differs from the Galerkinprocedure in that it is formulated from an implicit, discrete-time representation of Eq. 1, as opposedto the continuous-time representation for the Galerkin method [32]. For exemplary purposes, weconsider the first-order backward Euler scheme for the temporal discretization. Extension to otherimplicit methods is straightforward.The LSPG method is based upon the fully discretized algebraic version of Eq. 1, Q np − Q n − p ∆ t = ˆ R (cid:0) Q np , Ω nref (cid:1) , (7)where n represents the physical time step. Following the same procedure as Eq. 5, Eq. 7 is projectedonto the k th test function, W nk , (cid:42) W nk , Q np − Q n − p ∆ t − ˆ R (cid:0) Q np , Ω ref ( t ) (cid:1)(cid:43) P = r n . (8)7here the right-hand-side term, r n , represents the residual of the reduced system equations. Insteadof using the same POD basis for the test function as in the Galerkin method, Eq. 5, the goal ofthe LSPG method is to find a test function, W nk , such that the right-hand-side residual in Eq. 8, r n , is minimized. Similar to Least squares finite element methods [35], this minimization providesnumerical stabilization.By solving the minimization problem, a new time-varying test function, W nk , can be defined as, W nk = P − Φ k ∆ t − (cid:32) ∂ ˆ R ∂ Q p (cid:33) P − Φ k , (9)where (cid:16) ∂ ˆ R ∂ Q p (cid:17) is evaluated using the POD expansion, Eq. ?? [32]. Substituting Eq. 9 into Eq. 8, thediscretized ODE system can be obtained for the LSPG method, ˆ M n a n − a n − ∆ t = ˆ f (cid:0) a n , Ω nref (cid:1) , (10)where the elements in matrix in ˆ M n are, ˆ m k,j = (cid:104) W nk , Φ j (cid:105) P and those in ˆ f are, ˆ f k = (cid:68) W nk , ˆ R (cid:16) Q np , Ω nref (cid:17)(cid:69) P with Q p approximated using Eq. ?? .The symmetric nature of the matrix, ˆ M n , improves robustness from a global stability standpoint.The LSPG method requires implicit time integration, as indicated in Eqs. 7 and 10, while thestandard Galerkin method can be solved using either explicit or implicit schemes. Applying theLSPG method to Eq. 1 with explicit time discretization results in a ROM ODE system identical tothe Galerkin result. In the section below, we assess both Galerkin and LSPG methods in terms ofROM stability. III. Test Problem for Reduced Order Model Assessment
To assess the capabilities of ROMs for representing realistic combustion flowfields, a 2D-planarrepresentation of a generic laboratory-scale combustor [36] is used. This simplified model allowsROM capabilities to be evaluated while maintaining the essential physics of interest. The specificconfiguration is shown in Fig. 1. The problem consists of a shear coaxial injector with an outerpassage, T , that introduces fuel near the downstream end of the inner passage and a coaxial centerpassage, T , that feeds oxidizer to the combustion chamber. The T stream contains gaseous methanewhile the T2 stream is 42% gaseous O and 58% gaseous H O .8 ig. 1 Overview of 2D planar rocket injector benchmark problem. Operating conditions are maintained similar to conditions in the laboratory combustor [36, 37]with an adiabatic flame temperature of approximately 2700K and an imposed chamber pressure of1.1MPa. Both the T and T streams are fed with constant mass flow rates; 5.0kg/s and 0.37kg/srespectively. A non-reflective boundary condition is imposed at the downstream end to controlacoustic effects on the combustion dynamics. For all the cases in the current paper, a 10% sinusoidalperturbation at 5000Hz is imposed at the downstream boundary to generate FOM solutions forPOD mode generation. Combustion is represented by the single-step global model of Westbrookand Dryer [38], CH + 2 O → CO +2 H O. (11)As reported in Ref [18], stable, accurate reconstruction of CFD solutions of flows with stiff chemistryis highly challenging. To diminish these difficulties, the present simulations are based upon a reducedpre-exponential factor, A = 2 × , which is a factor of ten smaller than the value in [38] andcorresponds to a characteristic chemical time scale of approximately 0.8 µs and a laminar flamethickness of approximiately 1 mm . The FOM was calculated using a constant time step at 0.1 µs andgrid size, 0.18 mm , in the reacting regions. Even with this reduced reaction rate, the resulting ROMsremain highly temperamental and provide a clear example of the additional difficulties engenderedwhen reactions are present. This reduced reaction rate provides conditions more favorable for ROMdevelopment than the original stiff value while maintaining representative flame dynamics in thecombustor.A representative instantaneous snapshot of the reacting FOM solutions is shown in Fig. 2 todemonstrate the overall character of the flowfield and to highlight the dominant physics in theproblem of interest. The combustion dynamics are characterized by highly dispersed pockets of9 ig. 2 Representative instantaneous snapshots of temperature and heat release from CFDsimulation. intense heat release that are intermittently distributed in both space and time. The temperatureand heat-release contours span a wide range of scales from the small eddies in the shear layers to thelarge-scale recirculation zone behind the dump plane. All these unique features introduce varyinglevels of difficulty in constructing a robust ROM. IV. Numerical Investigations
As indicated above, ROMs for reacting flows tend to be less robust than those for non-reactingflows. The goal of the paper is to identify the main underlying causes for the challenges arising indeveloping ROMs for reacting flows. Having identified these issues, one can then use this under-standing to improve ROMs for combustion dynamics applications and potentially for non-reactingapplications as well. To this end, we conduct detailed investigations on the test problem illustratedin Fig. 1. To elucidate the source of ROM stability issues for reacting flows, we consider bothreacting and non-reacting versions of the test problem.The FOM time step for all calculations was 0.1 µs , and data were stored over a time duration of1.0ms, corresponding to a total of 10,000 time steps. The ROM database for POD mode generationwas established by down-sampling the FOM solution every 50 time steps resulting in a total of200 POD snapshots. Clearly, both the length of the time sample and the sampling interval areimportant issues in any reduced-order representation. The duration of the time sample dictatesthe lowest frequency that can be represented by the ROM, while the sampling interval restricts the10igh-frequency content. The above values were chosen as a compromise to represent the dominantphysics in the present test problem. A sensitivity study on the sampling interval is presented inAppendix B. A. Energy Content as a Function of the Number of POD Modes
The POD characteristics of the reacting flow problem are first investigated to understand howwell the POD modes represent the original CFD dataset. The representation is based upon thePOD residual energy:POD Residual Energy ( N p ) , %= (cid:32) − (cid:80) N p n =1 σ n (cid:80) N p,total n =1 σ n (cid:33) × , (12)where N p is the number of POD modes included, and N p,total (= 200) is the total number of snapshotsin the dataset . The residual energy as a function of N p shown Fig. 3 reveals the amount ofinformation omitted by the POD representation for any particular number of modes. The resultsshow that approximately the first 15 modes must be included to capture 90% of the total energywhile 45 modes recover approximately 99% and at least 133 modes (66.5% of the complete set) areneeded to retrieve 99.99% of the total energy. ROM applications in the literature are often presentedfor problems exhibiting less physical complexity such that they are able to capture 99.9% of theenergy with a dozen or so modes. The relatively large number of modes needed here is a testimonyto the wide range of dominant scales in the present problem as noted in Fig. 2. B. ROM Reconstruction Errors as a Function of the Number of POD Modes
The above results illustrate the number of modes required to represent a given fraction of theenergy in the original CFD solution. In the present section, we examine the accuracy with whichROM solutions can reproduce the original signal. Results are presented for ROMs based on bothGalerkin and least squares Petrov-Galerkin (LSPG) projection. Galerkin schemes are commonlyemployed in reduced-order modeling, but it is well-known that they are susceptible to instability innon-symmetric systems, especially when details are under-resolved. The symmetrized ODE systemprovided by the LSPG method should provide stability improvements and more robust ROMs. Suchimprovements have been demonstrated for non-reacting flows [32]. Here, we provide a three-way11 ig. 3 POD residual energy distribution from 2D rocket injector simulation. comparison of the methods: the Galerkin procedure evaluated for both explicit (RK-4) and implicittime marching of the ODE system, and the LSPG system with implicit integration. (As notedabove, with explicit integration the LSPG and Galerkin systems become identical.)Comparisons of the three methods are summarized in Fig. 4 in terms of the ROM reconstructionerror, Global ROM Reconstruction Error ( n p ) = 1 n nvar n var (cid:88) j =1 (cid:32) N I NI (cid:88) i =1 (cid:13)(cid:13)(cid:13)(cid:13)(cid:13) − q (cid:48) rom,j,i ( t ) q (cid:48) pod,j,i ( t ) (cid:13)(cid:13)(cid:13)(cid:13)(cid:13) (cid:33) (13)where q (cid:48) pod,j,i ( t ) = P − (cid:80) n p k =1 a pod ( t ) φ k,j,i represent the particular solution variable. Note the errorsin all n var variables are computed independently, summed and averaged. To ensure time step doesnot have a major effect, all three ROM versions are temporally integrated with three differenttime steps, 0.2, 0.1 and 0.05 µs . The 0.1 µs value coincides with the time step used in the FOMcalculation.Figure 4 presents the reconstruction error for all three ROM methods with results plotted formode number increments of 20 POD modes. Using more than 120 modes consistently results instable ROMs for all three methods whereas less than 100 modes lead to unstable ROMs. Using thereconstruction error as the sole basis for evaluation, there appears to be little reason for preferringone method over the other. In particular, the expected improvement with the LSPG method is12 ig. 4 Global ROM reconstruction error comparisons between different projection-based ROMmethods with varied computing time steps for the 2D benchmark reacting flow simulation. not realized in the present reacting flow calculations. As indicated later, this is probably becausethe major challenges in reacting flow ROMs arise from under-resolved local small-scale dynamicscorresponding to flame dispersion and high temperature/species gradients and these completelyoverwhelm the global, system-level stabilization provided by the LSPG method [19]. Later calcula-tions in which the local errors are controlled show the LSPG method produces a slight improvementover the Galerkin method.The results in Fig. 4 also indicate that all three methods are quite insensitive to the time step.Changing the time-step has a very minor impact on the LSPG solutions and an even smaller effecton the Galerkin predictions. Carlberg et al. [32] have identified a potential sensitivity to time stepin the LSPG method, but it does not appear to be a significant issue in the present work. Settingthe ROM time step equal to the FOM time step gives essentially identical reconstruction errorsfor all three methods. The minor differences in the LSPG method show that the ROM solutionswith the larger time step (0.2 µs ) are marginally less accurate, but that the (0.05 µs ) and (0.1 µs ) solutions are essentially identical. An issue that may be a factor in the present results isthat implicit time stepping (for either the Galerkin or LSPG methods) requires the calculation ofJacobians from the ROM-reconstructed solutions. It is possible that accuracy may be impacted notonly by the truncation of POD modes (i.e., small-scale dynamics) but also by errors in the computed13 ig. 5 Global ROM reconstruction error comparisons between different projection-based ROMmethods with computing time step, ∆ t = 0.1 µs , for the 2D benchmark reacting flow simula-tion. Jacobians. In general, Jacobian calculations are more challenging for reacting flow simulations withstiff kinetics, and these errors may escalate in regions where the solution strays from the originalcalculation. This issue is highlighted later as a potential issue in the implicit Galerkin results. Theoverall insensitivity to time step indicates the capability of ROMs in providing consistent modelingfidelity, even when larger time steps are used, thereby resulting in improved computational efficiency.To investigate the effect of the number of modes on stable ROMs in more detail, Fig. 5 showsreconstruction results for mode numbers between 100 and 120 POD modes based upon a modenumber increment of two. For consistency with the FOM, the time step in all ROM calculationshas been selected as ∆ t = 0.1 µs . As highlighted in the zoomed-in view in Fig. 5, at least 104 PODmodes are required to generate a stable ROM for the Galerkin method with explicit time marchingand at least 110 modes for both the LSPG and Galerkin methods with implicit time marching.More importantly, using two fewer POD modes for any of the three methods leads to unstableROMs indicating a strong sensitivity to the number of POD modes.These observations naturally lead to the following questions:1. Why does a difference of two POD modes have a significant impact on ROM stability?14. What are the underlying factors that determine ROM stabilization in reacting flow simula-tions?To address these issues, more detailed investigations are presented in the following sections. Inparticular they focus on three potential factors that may impact numerical stability of the ROMs:1) conservation law preservation; 2) loss of dissipation, and 3) the existence of unphysical localphenomena. C. Conservation Law Preservation
Carlberg et al. [39] have pointed out that neither Galerkin nor LSPG projection methods im-plicitly preserve fluid conservation laws for finite-volume models and have demonstrated that ROMstability can be improved by explicitly enforcing conservation in the computed ROMs. Therefore, asa first evaluation of potential contributing factors to ROM instabilities, the degree to which primaryand secondary conservation laws are violated in the ROM solutions is assessed. As an indicator ofprimary conservation law preservation, we define a global conservation error as:
Global Conservation Error ( n p ) = n eq (cid:80) n eq j =1 (cid:32) (cid:82) V ∂Qj∂t dV + (cid:80) NBCk =1 (cid:82) ∂ Ω k ( F − F v ) j · dS − (cid:82) V H j dV (cid:80) NInflowBCk =1 (cid:82) ∂ Ω k ( F − F v ) j · dS (cid:33) (14)where j represents the particular equation. In determining the conservation error, the temporalstorage, flux and source terms for each governing equation (continuity, momentum, energy andspecies) are computed for each cell at each time step and globally integrated to document theconservation errors. The errors are then normalized by the fluxes at the inflow boundaries andsummed over all governing equations. The temporal evolution of the conservation errors in theLSPG and the explicit and implicit Galerkin ROM solutions are then compared with each otherand the errors in the FOM. The results are summarized in Fig. 6.Conservation errors for each ROM method are presented for four different numbers of modes:the highest unstable and lowest stable number of modes for each method, plus 120 and 200 modesfor all three methods. As noted in Fig. 5, the lowest-unstable/highest-stable number of modes forthe Galerkin explicit method is 108 and 110, while the corresponding pair for the Galerkin implicitand LSPG implicit methods is 102 and 104. The 120-mode case provides a common assessment15 ig. 6 Global primary conservation error comparisons between different projection-basedROM methods with varied computing time steps for the reacting flow simulations. of all three methods at a number of modes where all methods perform well. The 200-mode casecorresponds to the inclusion of the entire POD database in the ROM.The FOM solutions exhibit conservation errors on the order of − (Fig. 6), a level consideredthat is typically considered sufficient for reacting flow simulations. The conservation errors in theROM solutions are more than an order of magnitude larger than the FOM errors, but in all cases,the error is decreased as additional POD modes are included. Even though the primary conservationerrors in the ROMs are larger than those in the FOM, the stability of the ROMs does not appear tobe related to an increase in conservation error, especially comparing the mode numbers straddlingthe unstable/stable boundary. For example, the 102-mode ROM from the explicit time marchingGalerkin solution exhibits a very similar error evolution as the 104-mode ROM up until the solutionblows up ( t ∼ s ). Similar trends are observed for the Galerkin implicit and LSPG methodsalthough the solution blows up at somewhat different times for each method. It is of particular notethat the errors remain nominally constant until blow-up occurs and then exhibit large increases.The preservation of secondary conservation is also of importance in numerical schemes [40] andits importance in ROM stability has been emphasized by Afkham et al. [41]. Both the FOM andROM solutions are based on the primary conservation laws, but secondary conservation (kinetic16 ig. 7 Global secondary conservation error comparisons between different projection-basedROM methods for the reacting flow simulations. energy, k , and entropy, s ) is evaluated by means of the following, Kinetic Energy:
DρkDt = u j ρu j Dt − k DρDt (15) Entropy:
DρsDt = 1 T (cid:32) D ( ρh − p ) Dt − DρkDt − (cid:88) k ( h k − T s k ) DρY k Dt (cid:33) (16)where k = u i . Specifically, kinetic energy ( k ) conservation is evaluated using mass and momentumconservation while entropy ( s ) conservation is obtained from energy, kinetic energy and speciesconservation. Similar to primary conservation errors, secondary conservation errors are averaged forkinetic energy and entropy for both FOM and ROMs and compared in Fig. 7. The FOM shows aone-order of magnitude larger error in secondary conservation than primary while the ROMs exhibitsimilar errors in both but still somewhat larger errors than in the FOM. Again, there appears tobe no direct connection between secondary conservation errors and ROM stability, suggesting thatconservation errors are not a dominant contributor to the stability issues observed in Figs. 4 and 5. D. Loss of Dissipation
A second potential contributor to ROM instability is that the amount of dissipation in the ROMODE system may be less than that in the FOM PDE because high-energy modes that contributeto dissipation are truncated [23]. Introducing additional dissipation via a closure model [23] or17umerical dissipation [7] has been demonstrated to improve ROM stability in non-reacting flowproblems. Accordingly, we compare the level of artificial dissipation (AD) in the ROM solutionswith that in the FOM.Fractional change of AD = (cid:82) V |∇ · (cid:126)F d | ROM − |∇ · (cid:126)F d | F OM dV (cid:82) V |∇ · (cid:126)F d | F OM dV × (17)where ∇ · (cid:126)F d = V (cid:80) N face n =1 (cid:104) Γ p | Γ − p A p | ( ˆ Q p,R − ˆ Q p,L ) (cid:105) · (cid:126)A n based upon a Roe-finite-volume, spatialdiscretization. A positive fractional change implies the AD in the ROM is larger than that in theFOM, while a negative value implies the converse. As with the conservation-error comparisons, thelevel of artificial dissipation in the ROMs is computed for different numbers of modes for the threemethods. The results are shown in Fig. 8.The first observation from Fig. 8 is that the relative change in AD is positive for all ROMs.For the present calculations, all ROMs contain more artificial dissipation than the FOM (As notedlater, the increased AD in the ROM solutions may be a result of the steeper gradients in the ROMsolutions compared to those in the FOM solutions). As additional POD modes are included, theoverproduction of AD in the ROMs gets smaller, but even the 200-mode ROM contains some 10- 15% more AD than the FOM, with the Galerkin explicit being slightly higher than the LSPGand Galerkin implicit methods. Second, for all three methods, stable ROMs are obtained even forAD overproduction levels up to 80%. More importantly, the level of AD in the unstable ROMs(the blue line in each plot) does not behave significantly different from that of the stable ROMs inwhich only two more POD modes are included. Though an excess accumulation of AD is present inthe Galerkin implicit scheme before it blows up, the ROMs of the other two methods blow up at atime where the AD is near a local minimum. The same assessment has been performed for viscousdissipation (not shown) and leads to similar conclusions.In summary, there appear to be no distinguishable connection between ROM blow-up and thelevel of AD or the level of primary or secondary conservation. A remaining possibility is that insta-bilities in the ROM solutions arise because the ROM solutions generate unphysical local phenomena.This issue is discussed in the following section. 18 ig. 8 Evaluation of artificial dissipation (AD) overproduction between different projection-based ROM methods for the 2D benchmark reacting flow simulation.E. Unphysical Local Phenomena Thus far, the ROM stability contributor assessment is based on evaluations of global properties,similar to most ROM stabilization methods in the literature including LSPG. In this section, thefocus of the assessment is on local phenomena. Due to the large number of degrees of freedom inthe FOM and ROM solutions (10,000 snapshots ×
40k cells), the minimum and maximum values ofthe primitive variables (pressure, velocities, temperature and species) are first investigated for eachsnapshot from the ROMs to help identify local unphysical phenomena. These investigations drawour attention to the temporal evolution of minimum temperature in the domain at each time stepas presented in Fig. 9 for the three different ROM methods. The minimum temperatures in theFOM solution are also included for reference.An immediate observation from Fig. 9, is that all three ROM methods contain temperaturesthat are far below the temperature of the incoming fluids whereas the minima in the FOM remainsat essentially the incoming temperature. The full 200-mode ROMs show reasonable temperaturelevels that correlate well with the FOM solution, but as the number of modes in any of the threeROM methods is decreased, the minimum temperature decreases to unphysical values that aresubstantially below the incoming fluid temperature. The minimum temperature fluctuates in timein response to the imposed periodic forcing at the downstream end, but decreases consistently as19 ig. 9 Global minimum temperature values for projection-based ROM methods for the react-ing flow simulation. time goes on. Most significantly, the unstable ROM for each of the three methods (102 modes forGalerkin with explicit time marching and 108 modes for the LSPG and Galerkin implicit methods)decrease to absolute zero (or below). This ultimately unphysical temperature is clearly the causefor the abrupt solution departure. Unphysical temperatures prove to be the major reason unstableROMs are generated in this reacting flow problem.Even before negative temperatures are reached, the presence of unphysically low temperaturescan cause other problems in the ROM calculations. Unphysical temperatures lead to errors indensity, sound speed and thermal properties, which can give rise to ill-conditioned Jacobian matricesin the ROM calculations causing further errors. Of particular importance in reacting flow problems isthe strong interaction between the local temperature and the highly nonlinear Arrhenius terms thatappear in the source terms in the species equations. A small error in temperature generates a muchlarger error in the local reaction rate, thus exacerbating the error. The presence of unphysical localtemperatures addresses both of the questions proposed at the end of Section IV B and explains whyLSPG does not provide significant stability improvements in the resulting ROMs. The unphysicaltemperatures appear in both Galerkin and LSPG solutions so the global stability improvements ofthe LSPG method are swamped by these local errors.To provide more insight into the source of these unphysical temperatures, instantaneous span-20ise temperature profiles for the FOM and ROM solutions are compared at the approximate stream-wise location where the minimum temperature is observed. The temperature excursions for all threeROM methods are similar but for brevity, only results from Galerkin ROMs with implicit time step-ping are shown in Fig. 10. The axial location shown is x = 0.009m (see Fig. 1) and the time is chosento coincide with the time step at which the 108-mode ROM solution in Fig. 9 ( ∼ F. Temperature Fluctuations in Non-reacting Flows
As a means of distinguishing the effects of steep temperature gradients from the effects of chem-ical reactions, we briefly consider a set of non-reacting flows. For consistency with the reacting flowproblem, three sets of FOM simulations of non-reacting flow are performed using the configuration21 ig. 10 Spanwise temperature profile comparisons between FOM and ROMs at representativestreamwise location (x = 0.009m in Fig. 1) and blow-up time instance (identified from Fig. 9)for Galerkin ROM with implicit scheme. in Fig. 1. The incoming fluids contain the same species concentrations in the T and T streams asthe reacting flow case but three different temperature ratios are considered in the two streams. Inall cases, the fluid in stream T is fixed at 300K but the temperature ratios are chosen as T / T = 1.0, 3.0 and 4.0) to investigate the effects of temperature variations. Combustion is omitted bysetting the pre-exponential factor to zero.Instantaneous snapshots from the three FOM solutions are compared in Fig. 11. Overall, therepresentative dynamics of all cases are qualitatively similar to the reacting flow case and showmany of the characteristics expected in reacting flow. Each flowfield exhibits small-scale eddiesoriginating within the recess between the exit of the T1 stream and the dump plane. These smalleddies are then shed into the combustor and follow the large-scale recirculation zone created by thedump plane, eventually being absorbed into the larger scales.The residual energy in the three POD datasets is shown as a function of the number of modesin Fig. 12. Comparison with Fig. 3 reveals that the amount of information omitted by the POD22 ig. 11 Representative instantaneous snapshots comparisons for 2D benchmark non-reactingflow simulations with varied temperature ratios.Fig. 12 POD residual energy distribution comparisons for 2D benchmark non-reacting flowsimulations with varied temperature ratios. representation for any particular number of modes is essentially the same as in the reacting flow case,independent of the temperature ratio in the non-reacting flow solution with only modest differencesbetween the four datasets.ROM reconstruction errors are plotted as a function of the number of POD modes in Fig. 13for all three temperature ratios. For the T / T =1.0 case, continuous monotone convergence from20 to 200 modes is observed and all generated ROMs are stable. The errors in the two non-unity23 ig. 13 Global ROM reconstruction error comparisons for the 2D benchmark non-reactingflow simulation with varied temperature ratios (unstable ROMs marked in hollow symbols). temperature ratio cases, however, are very large when a small number of modes is included but dropabruptly after 60 modes in the T / T = 3.0 case and after 80 modes in the T / T = 4.0 case.ROMs constructed with less than these numbers of modes are unstable for both cases. Increasing thetemperature ratio requires more POD modes to reach a stable ROM and makes it more challengingto generate robust ROMs. Comparison with Fig. 4 shows that the T / T = 3.0 and 4.0 cases arehighly analogous to the reacting flow cases. Inspection of these non-reacting flows again shows strongtemperature over-and under-shoots near steep temperature gradients. The complexities generatedby steep temperature gradients is one of main challenges in ROM development for reacting flowsimulations. G. Implementation of Temperature Constraints
Both the reacting and non-reacting solutions suggest that the major cause of ROM stabilityarises because of Gibbs phenomena near steep temperature gradients leading to unphysical temper-atures that are not present in the CFD database. As a verification that the temperature excursionsare the source of instability, simple temperature constraints are imposed on the ROM calculations.After each ROM time iteration, the temperatures are constrained in between 250 and 2750K. The24ower bound of the constraints, 250K, is determined by the cold reactant temperature, 50K belowthe incoming temperature of 300K while the higher bound, 2750K, is 50K above the adiabatic tem-perature of the flame of 2700K. Temperature constraints are also imposed for ROMs of the twohigh temperature ratio non-reacting flow cases in Fig. 13 with the same lower bound but 950K and1250K for the higher bounds for T / T = 3.0 and T / T = 4.0 respectively. (Adding constraintson the T / T = 1.0 case has no impact as the limits are never exceeded.)The global ROM reconstruction errors with temperature constraints imposed are compared inFig. 14 with solutions without the constraints. The plots at the top compare the non-reacting cases.Results are shown for T / T = 3.0 and 4.0 with and without constraints imposed, and for T / T = 1.0 (where constraints have no effect). Imposing constraints brings the convergence error forboth T / T = 3.0 and 4.0 nearly into coincidence with that for the T / T = 1.0 case. The ROMconvergence error starts at a stable level and decreases monotonically as additions modes are added.All ROMs generated from 20 to 200 modes are stable and robust.The plots at the bottom of Fig. 14 present ROM convergence for the reacting flow problem forthe Galerkin and LSPG methods with and without temperature constraints. The Galerkin explicitand LSPG methods show substantial improvements when temperature constraints are employedbut the Galerkin implicit method shows little difference. The effect of temperature constraints ismost beneficial for the LSPG results. With temperature constraints, the LSPG method producesstable ROMs for all POD mode numbers from 20 to 200 with nearly monotonic convergence in ROMerrors. Without constraints, 110 modes are required for stability. Both Galerkin solutions retainregions of instability when small numbers of POD modes are used. For the Galerkin solution withexplicit time marching, stable ROMs can be obtained by including 60 or more POD modes whereaswithout temperature constraints a minimum of 104 modes is required. The Galerkin solutions withimplicit time marching, however, show essentially no improvement with temperature constraints.Implicit Galerkin solutions still require 110 POD modes to obtain stable ROMs. Further analysisof the Galerkin implicit results indicate that even though the existence of extreme local phenomenahas been prevented, the mass matrix of the ROM ODEs (similar to ˆ M n in Eq. 10) still exhibitspoor numerical properties that lead to numerical divergence during ROM iterations and eventually25 a) Non-reacting Flow(b) Reacting Flow Fig. 14 Global ROM reconstruction error comparisons for the 2D benchmark non-reactingand reacting flow simulations with temperature constraints (unstable ROMs marked in hollowsymbols). terminate the calculations. On the other hand, both Galerkin with explicit scheme and LSPGmethods produce symmetric mass matrices (identity matrices for the explicit scheme), which ingeneral have better numerical properties.As a final observation, a further benefit of temperature constraints is that they enable the ROMsto be extended to future-state predictions beyond the time interval used for basis construction.Comparisons between the FOM and 200-mode ROMs for reacting flow with and without temperature26 ig. 15 Evaluations of future-state predictive capabilities for LSPG ROM with and withouttemperature constraints based on pressure signal taken near the dump plane (x = 0 m in Fig.1). constraints are shown in Fig. 15 for the LSPG ROM methods. Future state predictions with theGalerkin methods are similar. Both with and without temperature constraints, the 200 POD modesproduce stable and accurate ROMs within the reconstruction interval (0.015 to 0.016s), however,as demonstrated in the pressure signal comparisons in Fig. 15, only the ROMs with temperatureconstraints are able to predict solutions beyond 0.016s. ROMs without constraints blow up shortlyafter the calculations proceed outside the reconstruction interval. Although the capability of future-state predictions requires further investigations and improvement, as is the case in simpler non-reacting flow problems [32]. It should also be pointed out that such advancement can also beachieved for Galerkin method as well, both explicit and implicit with the temperature constraints.Similar improvements in future-state predictions were achieved for the non-reacting flow cases.
V. Conclusion
The accuracy, stability and robustness of reduced order models (ROM) for reacting flow ap-plications has been investigated in a detailed manner. The goal is to reveal the prime underlyingcauses for numerical instability and evaluate corresponding resolutions to improve the robustnessof the ROMs. Galerkin and least squares Petrov-Galerkin (LSPG) methods were compared for27 representative combustor problem to assess the baseline ROM characteristics in reacting flowsimulations. Overall, ROM stability appears to be highly sensitive to the number of POD modesincluded and the addition of only two modes can switch from unstable to stable ROMs. Attentionwas focused on investigating three potential contributing factors for ROM stability: 1) conservationlaw preservation; 2) loss of dissipation; and 3) unphysical local phenomena.In addressing the first potential reason for ROM sensitivity, preservation of primary and sec-ondary conservation was quantified and investigated for FOM and ROMs. ROMs exhibit conser-vation errors more than an order of magnitude higher than those in the FOM, but they remainsufficiently small that stable ROMs are still be produced. Further, the level of conservation errorsin unstable ROMs is similar to that in stable ones suggesting that the violation of conservation isnot a significant contributor to ROM stability.Next, the level of artificial dissipation in the ROMs was shown to be consistently larger thanin the FOM with increases of 10% to 80%. Again, the global dissipation level was observed to besimilar in stable and unstable ROMs, and does not appear to be a significant cause of ROM stability.Third, attention was directed toward unphysical local phenomena, which strongly suggest thatthe blow-up of ROMs is driven by the appearance of unphysical temperatures arising from Gibbsphenomena near steep temperature gradients. Similar conclusions were observed for non-reactingflows with steep temperature gradients.Based on the identification of unphysical temperatures as a major contributor to instability, theROM temperature fields were clipped to preclude unphysical values. With temperature constraints,both Galerkin and LSPG methods exhibited significant improvements in ROM stability and errorconvergence. Imposition of temperature constraints was also shown to enable effective future-stateprediction. It should be pointed out that the temperature constraints are not adhoc. The lowerbound is determined based on the cold reactant temperature while the higher bound is set based onthe adiabatic flame temperature. Clearly, more sophisticated limiters are preferred, but the presentresults highlight the necessity of controlling sharp gradients.Although the imposition of temperature constraints yields significant improvement in ROMperformance, it does not eliminate the Gibbs phenomena. This is indicative of the deficiency of flux28imiters in eliminating Gibbs phenomena in ROMs or the lack of appropriate closure or stabilizationto account for the truncated bases. Therefore, more comprehensive improvements can be achievedby developing ROM-specific limiters [42], more systematic closure models [43] and by using adaptivebases [18, 44] to eliminate Gibbs phenomena. Further consideration should be given to the challengesinvolved in improving the efficiency of the ROM using sparse sampling techniques [18, 45].
APPENDIX A: GOVERNING EQUATIONS FOR FULL ORDER MODEL
The full order model computations are carried out with an in-house CFD code, the GeneralEquations and Mesh Solver (GEMS), the capabilities of which has been successfully demonstratedin modeling rocket combustion instabilities [1]. GEMS solves the conservation equations for mass,momentum, energy and species mass fractions in a coupled fashion, ∂Q∂t + ∇ · (cid:16) (cid:126)F − (cid:126)F v (cid:17) = H, (A1)where Q is the vector of conserved variables defined as, Q = (cid:18) ρ ρu ρv ρw ρh − p ρY l (cid:19) T with ρ representing density, u , v and w representing velocity field, Y l representing the l th species massfraction and the total enthalpy h is defined as, h = h + ( u i ) = (cid:80) l h l Y l + ( u i ) .The fluxes have been separated into inviscid, (cid:126)F = F i (cid:126)i + F j (cid:126)j + F k (cid:126)k and viscous terms, (cid:126)F v = F v,i (cid:126)i + F v,j (cid:126)j + F v,k (cid:126)k . And the three inviscid fluxes are, F i = ρuρu + pρuvρuwρuh ρuY l , F j = ρvρuvρv + pρvwρvh ρvY l and F k = ρwρuwρvw + pρw + pρwh ρwY l (A2)29he viscous fluxes are, F v,i = τ ii τ ji τ ki uτ ii + vτ ji + wτ ki − q i ρD l ∂Y l ∂x , F v,j = τ ij τ jj τ kj uτ ij + vτ jj + wτ kj − q j ρD l ∂Y l ∂y and F v,k = τ ik τ jk τ kk uτ ik + vτ jk + wτ kk − q k ρD l ∂Y l ∂z (A3)where D l is defined to be the diffusion of the l th species into the mixture. In practice, this is anapproximation used to model the multicomponent diffusion as the binary diffusion of each speciesinto a mixture.The heat flux in the i th direction, q i , is defined as, q i = − K ∂T∂x i + ρ N (cid:88) l =1 D l ∂Y l ∂x i h l + Q source (A4)The three terms in the heat flux represent the heat transfer due to the conduction, species diffu-sion and heat generation from a volumetric source (e.g. heat radiation or external heat source)respectively.The shear stress, τ , is also found in the viscous flux and defined in terms of the molecularviscosity and velocity field, τ ij = µ (cid:18) ∂u i ∂x j + ∂u j ∂x i − ∂u m ∂x m δ ij (cid:19) (A5)The source term, H includes a single entry for each of the species equations signifying theproduction or destruction of the l th species, ˙ ω l , which is determined by the chemical kinetics [38], H = (cid:18) ω l (cid:19) T (A6) APPENDIX B: INVESTIGATION ON THE EFFECTS OF SNAPSHOT SELECTION
With the improvement on ROM stability via the imposition of temperature constraints asdemonstrated in Section IV G, further sensitivity studies are carried out to investigate the impactof snapshots selection on ROM characteristics. Accordingly, the FOM solutions have been down-sampled every 50, 25 and 16 time steps for POD mode generation, which results in 200, 400 and626 snapshots respectively. It should be noted that all the studies above were carried out using a30own-sampling rate of 50 time steps, and by increasing the sampling rates, it is expected that morehigh-frequency dynamics will be included in the ROM evolution. The comparisons of the PODresidual energy (Eq. 12) and the global ROM construction error (Eq. 13) are shown in Fig. 16.Fig. 16a reveals the amount of information omitted by the POD representation for any particularnumber of modes. It is readily seen that the residual energies for all three sampling rates follownearly identical decay rates for approximately the first 60 POD modes. Concurrently, it is noted thatthese first 60 modes contain nearly 90% of the total energy. The relatively large fraction of modesrequired to recover 90% of the energy, in contrast to non-reacting flow problems in the literature,emphasizes the wide range of important scales in the present problem. When more than 60 modesare included in the POD data set, substantial differences appear between the three sampling rates asmore high-frequency dynamics (possibly including numerical noise in FOM solutions) are includedwith the increase of sampling rates, which usually accumulates in the low-energy POD modes.The ROMs are constructed using LSPG with temperature constraints introduced in Sec-tion IV G, and the reconstruction errors are compared in Fig. 16b for three sampling rates. Itcan be readily seen that the ROMs from all three sampling rates are stable and exhibit largelymonotone error convergence. More importantly, as the sampling rate increases, the accuracy of theROM reconstruction appears to be consistently improved by at least a factor of 2 to 3, comparingthe 200-mode ROMs between sampling rate = 50 and 25 and the 400-mode ROMs between samplingrate = 25 and 16. This is expected since with a higher sampling rate, there is less information loss,and therefore, the resulting ROMs are anticipated to provide a more accurate representation of theFOM solutions.
ACKNOWLEDGMENTS
The authors acknowledge support from the Air Force under the Center of Excellence grantFA9550-17-1-0195, titled Multi-Fidelity Modeling of Rocket Combustor Dynamics. The authorsalso would like to acknoledge the computing resources provided by the NSF via grant 1531752 MRI:Acquisition of Conflux, A Novel Platform for Data-Driven Computational Physics (Tech. Monitor:Stefan Robila). 31 a) POD residual energy distribution (b) Global ROM reconstruction error (LSPG withtemperature constraints)
Fig. 16 Parametric investigations on the effects of down-sampling rates.REFERENCES [1] Harvazinski, M. E., Huang, C., Sankaran, V., Feldman, T. W., Anderson, W. E., Merkle, C. L., andTalley, D. G., “Coupling between hydrodynamics, acoustics, and heat release in a selfexcited unstablecombustor,”
Physics of Fluids , Vol. 27, 2015, pp. 1–27.[2] Staffelbach, G., Gicquel, L. Y. M., Boudier, G., and Poinsot, T., “Large Eddy Simulation of self excitedazimuthal modes in annular combustors,”
Proceedings of the Combustion Institute , Vol. 32, No. 2, 2009,pp. 2909–2916,doi:http://dx.doi.org/10.1016/j.proci.2008.05.033.[3] Wolf, P., Balakrishnan, R., Staffelbach, G., Gicquel, L. Y. M., and Poinsot, T., “Using LES to StudyReacting Flows and Instabilities in Annular Combustion Chambers,”
Flow, Turbulence and Combustion ,Vol. 88, No. 1, 2012, pp. 191–206,doi:10.1007/s10494-011-9367-7.[4] Gicquel, S. G., L.Y.M. and Poinsot, T., “Large Eddy Simulations of gaseous flames in gas turbinecombustion chambers,”
Progress in Energy and Combustion Science , Vol. 38, No. 782-817.[5] Lumley, J. L. and Poje, A., “Low-dimensional models for flows with density fluctuations,”
Physics ofFluids , Vol. 9, No. 7, 1997, p. 2023,doi:10.1063/1.869321.[6] Graham, W., Peraire, J., and Tang, K., “Optimal Control of Vortex Shedding Using Low Order Models art I: Open-Loop Model Development,” International Journal for Numerical Methods .[7] Lucia, D. J. and Beran, P. S., “Projection methods for reduced order models of compressible flows,”
Journal of Computational Physics , Vol. 188, No. 1, 2003, pp. 252–280,doi:10.1016/s0021-9991(03)00166-9.[8] Taira, K., Hemati, M. S., Brunton, S. L., Sun, Y., Duraisamy, K., Bagheri, S., Dawson, S., and Yeh,C.-A., “Modal analysis of fluid flows: Applications and outlook,” arXiv preprint arXiv:1903.05750 .[9] Barbagallo, A., Dergham, G., Sipp, D., Schmid, P. J., and Robinet, J.-C., “Closed-loop control ofunsteadiness over a rounded backward-facing step,”
Journal of Fluid Mechanics , Vol. 703, 2012, pp.326–362,doi:10.1017/jfm.2012.223.[10] Barbagallo, A., Sipp, D., and Schmid, P. J., “Closed-loop control of an open cavity flow using reduced-order models,”
Journal of Fluid Mechanics , Vol. 641, 2009, p. 1,doi:10.1017/s0022112009991418.[11] Barbagallo, A., Sipp, D., and Schmid, P. J., “InputâĂŞoutput measures for model reduction and closed-loop control: application to global modes,”
Journal of Fluid Mechanics , Vol. 685, 2011, pp. 23–53,doi:10.1017/jfm.2011.271.[12] Lieu, T. and Farhat, C., “Adaptation of Aeroelastic Reduced-Order Models and Application to an F-16Configuration,”
AIAA Journal , Vol. 45, No. 6, 2007, pp. 1244–1257,doi:10.2514/1.24512.[13] Lucia, D. J., Beran, P. S., and Silva, W. A., “Reduced-order modeling: new approaches for computa-tional physics,”
Progress in Aerospace Sciences , Vol. 40, No. 1-2, 2004, pp. 51–117,doi:10.1016/j.paerosci.2003.12.001.[14] Munipalli, R., Zhu, X., Menon, S., and Hesthaven, J., “Model Reduction Opportunities in DetailedSimulations of Combustion Dynamics,” , 2014,doi:10.2514/6.2014-0820.[15] Huang, X. and Baumann, W. T., “Reduced-Order Modeling of Dynamic Heat Release for Thermoa-coustic Instability Prediction,”
Combustion Science and Technology , Vol. 179, No. 3, 2007, pp. 617–636,doi:10.1080/00102200600969979.[16] Huang, C., Anderson, W. E., Merkle, C., and Sankaran, V., “Exploration of POD-Galerkin Methodin Developing a Flame Model for Combustion Instability Problems,” in “7th AIAA Theoretical FluidMechanics Conference,” , 2014,doi:10.2514/6.2014-3222.
17] Huang, C., Anderson, W. E., Merkle, C., and Sankaran, V., “Exploration of POD-Galerkin Techniquesfor Developing Reduced Order Models of the Euler Equations,” in “51st AIAA/SAE/ASEE Joint Propul-sion Conference,” , 2015,doi:10.2514/6.2015-4163.[18] Huang, C., Xu, J., Duraisamy, K., and Merkle, C., “Exploration of Reduced-Order Models for RocketCombustion Applications,” in “AIAA SCITECH Forum 2018,” , 2018.[19] Huang, C., Duraisamy, K., and Merkle, C., “Challenges in Reduced Order Modeling of Reacting Flows,”in “AIAA Propulsion and Energy Forum 2018,” , 2018.[20] Huang, C., Anderson, W., Merkle, C., and Sankaran, V., “Multifidelity Framework for Modeling Com-bustion Dynamics,”
AIAA Journal .[21] Xu, J., Huang, C., and Duraisamy, K., “Reduced-Order Modeling Framework for Combustor InstabilitiesUsing Truncated Domain Training,”
AIAA Journal , pp. 1–15,doi:10.2514/1.j057959.[22] Rempfer, D., “On Low-Dimensional Galerkin Models for Fluid Flow,”
Theoretical and ComputationalFluid Dynamics , Vol. 14, 2000, p. 75âĂŞ88.[23] Bergmann, M., Bruneau, C. H., and Iollo, A., “Enablers for robust POD models,”
Journal of Compu-tational Physics , Vol. 228, No. 2, 2009, pp. 516–538,doi:10.1016/j.jcp.2008.09.024.[24] Noack, B. R., Papas, P., and Monkewitz, P. A., “The need for a pressure-term representation in empiricalGalerkin models of incompressible shear flows,”
Journal of Fluid Mechanics , Vol. 523, 2005, pp. 339–365,doi:10.1017/s0022112004002149.[25] Rowley, C. W., “Model Reduction for Fluids using Balanced Proper Orthogonal Decomposition,”
In-ternational Journal of Bifurcation and Chaos .[26] Willcox, K. and Peraire, J., “Balanced Model Reduction via the Proper Orthogonal Decomposition,”
AIAA Journal , Vol. 40.[27] Rowley, C. W., Colonius, T., and Murray, R. M., “Model reduction for compressible flows using PODand Galerkin projection,”
Physica D: Nonlinear Phenomena , Vol. 189, No. 1-2, 2004, pp. 115–129,doi:10.1016/j.physd.2003.03.001.[28] Barone, M., Segalman, D., Thornquist, H., and Kalashnikova, I., “Galerkin Reduced Order Models forCompressible Flow with Structural Interaction,” , 2008,doi:10.2514/6.2008-612.
29] Barone, M., Kalashnikova, I., Segalman, D., and Heidi, K., “Stable Galerkin reduced order models forlinearized compressible flow,”
Journal of Computational Physics , Vol. 228, 2009, p. 1932âĂŞ1946,doi:10.1016/j.jcp.2008.11.015.[30] Barone, M. F., Kalashnikova, I., Brake, M. R., and Segalman, D. J., “Reduced Order Modeling ofFluid/Structure Interaction,” Report, Sandia National Laboratories, 2009.[31] Amsallem, D. and Farhat, C., “On the Stability of Projection-Based Linear Reduced-Order Models:Descriptor vs Non-Descriptor Forms,” in Quarteroni, A. and Rozza, G., eds., “Reduced Order Methodsfor Modeling and Computational Reduction,” Springer MS & A Series, Springer.[32] Carlberg, K., Barone, M., and Antil, H., “Galerkin v. least-squares Petrov-Galerkin projection in non-linear model reduction,”
Journal of Computational Physics , Vol. 330, 2017, pp. 693–734,doi:10.1016/j.jcp.2016.10.033.[33] Roe, P., “Approximate Riemann Solvers, Parameter and Difference Schemes,”
Journal of ComputationalPhysics , Vol. 43, 1981, pp. 357–372.[34] Barth, T. and Jespersen, D., “The design and application of upwind schemes on unstructured meshes,”in “27th Aerospace Sciences Meeting,” American Institute of Aeronautics and Astronautics, AerospaceSciences Meetings, 1989,doi:doi:10.2514/6.1989-36610.2514/6.1989-366. Doi:10.2514/6.1989-366.[35] Hughes, T. J., Franca, L. P., and Hulbert, G. M., “A new finite element formulation for computationalfluid dynamics: VIII. The Galerkin/least-squares method for advective-diffusive equations,”
Computermethods in applied mechanics and engineering , Vol. 73, No. 2, 1989, pp. 173–189.[36] Yu, Y., Sisco, J. C., Rosen, S., Madhav, A., and Anderson, W. E., “Spontaneous Longitudinal Combus-tion Instability in a Continuously-Variable Resonance Combustor,”
Journal of Propulsion and Power ,Vol. 28, No. 5, 2012, pp. 876–887,doi:10.2514/1.b34308.[37] Yu, Y.,
Experimental and Analytical Investigations of Longitudinal Combustion Instability in a Contin-uously Variable Combustor (CVRC) , Thesis, 2009.[38] Westbrook, C. K. and Dryer, F. L., “Simplified Reaction Mechanisms for the Oxidation of HydrocarbonFuels in Flames,”
Combustion Science and Technology , Vol. 27, No. 1-2, 2007, pp. 31–43,doi:10.1080/00102208108946970.[39] Carlberg, K., Choi, Y., and Sargsyan, S., “Conservative model reduction for finite-volume models,”
Journal of Computational Physics , Vol. 371,doi:10.1016/j.jcp.2018.05.019.
40] Tadmor, E., “Entropy stability theory for difference approximations of nonlinear conservation laws andrelated time-dependent problems,”
Acta Numerica ,doi:10.1017/S0962492902000156.[41] Afkham, B., Ripamonti, N., Qian, W., and Hesthaven, J., “Conservative Model Order Reduction forFluid Flow,” , 2018.[42] Ray, D. and Hesthaven, J., “Anartificialneuralnetworkasatroubled-cellindicator,”
JournalofComputa-tionalPhysics , pp. 166–191.[43] Parish, E., Wentland, C., and Duraisamy, K., “The Adjoint Petrov-Galerkin Method for Non-LinearModel Reduction,”
Preprint submitted to Journal of Computational Physics ,doi:10.2514/1.j057959.[44] Peherstorfer, B. and Willcox, K., “Online adaptive model reduction for nonlinear systems via low-rankupdates,”
SIAM Journal on Scientific Computing , Vol. 37, No. 4, 2015, pp. A2123–A2150.[45] Carlberg, K., Farhat, C., Cortial, J., and Amsallem, D., “The GNAT method for nonlinear modelreduction: effective implementation and application to computational fluid dynamics and turbulentflows,”
Journal of Computational Physics , Vol. 242, 2013, pp. 623–647., Vol. 242, 2013, pp. 623–647.