The conventional method to predict the onset of laminar-turbulent transition in convectively unstable boundary-layer flows is based on the logarithmic amplification ratio, the N-factor, of the linear instability waves. To calculate the N-factor, the flow variables are decomposed into a laminar basic state solution and the linear disturbances, which are assumed to be harmonic in time. The most commonly used linear stability analysis approaches include the locally parallel linear stability theory (LST) and the nonlocal, weakly nonparallel parabolized stability equations (PSE). However, these methods do not account for strong streamwise gradients that are encountered in several configurations of interest, such as those in the vicinity of roughness elements, steps, gaps, or corners. To compute the linear evolution of disturbances along such strongly nonparallel regions, the harmonic linearized Navier–Stokes equations (HLNSE) need to be solved. The discretization of the HLNSE for spanwise/azimuthally inhomogeneous laminar basic states yields a linear system of complex arithmetic with a leading dimension of the order of 107 to 108 even in relatively simple flows. A combined multithread and multiprocessor algorithm is implemented for the direct solution of such linear systems. Results for a supersonic boundary layer over a three-dimensional roughness patch show good agreement with experimental measurements when the evolution of the instability waves over the roughness patch is included via the HLNSE. Additionally, inflow-resolvent analysis based on the HLNSE for discrete-roughness-induced disturbances in the nose tip of a blunt cone at Mach 6 demonstrates the importance of including the disturbance amplification along the near vicinity of the roughness element and separation region.
A stencil-adaptive SBP-SAT finite difference scheme is shown to display superconvergent behavior. As proof of concept, applied to the linear advection equation, it has a convergence rate O(Δx4) in contrast to a conventional scheme, which converges at a rate O(Δx3).
The linear amplification of modal disturbances that lead to boundary-layer transition in two-dimensional/axisymmetric hypersonic configurations is strongly reduced by the presence of a blunt nose tip, and the mechanisms underlying the low Mack’s second-mode [Formula: see text]-factor values at the observed onset of transition over the cone frustum are currently unknown. Linear nonmodal analysis has shown that both planar and oblique traveling disturbances that peak within the entropy-layer experience appreciable energy amplification for moderate to large nose-tip bluntness. The present study extends the previous linear analysis by including the nonlinear effects. Specifically, the harmonic Navier–Stokes equations (HNSE) are solved with a fully implicit formulation and the Newton–Raphson method. The increased number of degrees of freedom for the nonlinear system presents difficulties for solution strategies based on direct solution of the linearized system. Such difficulties are overcome by using the generalized minimal residual method (GMRES) iterative method with a preconditioner corresponding to a simplified Jacobian without the cross-derivative terms. The HNSE solver is verified by comparing with nonlinear parabolized stability equation results for the nonlinear evolution of planar waves in an incompressible Blasius boundary layer and in a Mach 6 flow over a blunt cone. Finally, nonlinear nonmodal results are presented for planar traveling disturbances over a blunt cone configuration with reduced transition [Formula: see text] factor as measured in wind-tunnel experiments. The nonmodal analysis demonstrates that entropy-layer disturbances generated close to the nose tip can seed the amplification of higher-frequency Mack’s second-mode instabilities farther downstream and hence can lead to a reduction in the transition [Formula: see text] factor.
We present a new methodology for computing sensitivities in evolutionary systems using a model-driven low-rank approximation. To this end, we formulate a variational principle that seeks to minimize the distance between the time derivative of the reduced approximation and sensitivity dynamics. The first-order optimality condition of the variational principle leads to a system of closed-form evolution equations for an orthonormal basis and corresponding sensitivity coefficients. This approach allows for the computation of sensitivities with respect to a large number of parameters in an accurate and tractable manner by extracting correlations between different sensitivities on the fly. The presented method requires solving forward evolution equations, sidestepping the restrictions imposed by forward/backward workflow of adjoint sensitivities. For example, the presented method, unlike the adjoint equation, does not impose any I/O load and can be used in applications in which real time sensitivities are of interest. We demonstrate the utility of the method for three test cases: (1) computing sensitivity with respect to model parameters in the Rossler system (2) computing sensitivity with respect to an infinite-dimensional forcing parameter in the chaotic Kuramoto-Sivashinsky equation and (3) computing sensitivity with respect to reaction parameters for species transport in a turbulent reacting flow.
Provably stable flux reconstruction (FR) schemes are derived for partial differential equations cast in curvilinear coordinates. Specifically, energy stable flux reconstruction (ESFR) schemes are considered as they allow for design flexibility as well as stability proofs for the linear advection problem on affine elements. Additionally, the curvilinear metric split-form for a linear physical flux is examined as it enables the development of energy stability proofs. The first critical step proves, that in curvilinear coordinates, the discontinuous Galerkin (DG) conservative and non-conservative forms are inherently different-even under exact integration and analytically exact metric terms. This analysis demonstrates that the split form is essential to developing provably stable DG schemes on curvilinear coordinates and motivates the construction of metric dependent ESFR correction functions in each element. Furthermore, the provably stable FR schemes differ from schemes in the literature that only apply the ESFR correction functions to surface terms or on the conservative form, and instead incorporate the ESFR correction functions on the full split form of the equations. It is demonstrated that the scheme is divergent when the correction functions are only used for surface reconstruction in curvilinear coordinates. We numerically verify the stability claims for our proposed FR split forms and compare them to ESFR schemes in the literature. Lastly, the newly proposed provably stable FR schemes are shown to obtain optimal orders of convergence. The scheme loses the orders of accuracy at the equivalent correction parameter value c as that of the one-dimensional ESFR scheme. Crown Copyright (C) 2022 Published by Elsevier Inc. All rights reserved.
The entropy conservative, curvilinear, nonconforming, p-refinement algorithm for hyperbolic conservation laws of Del Rey Fernandez et al. (2019), is extended from the compressible Euler equations to the compressible Navier-Stokes equations. A simple and flexible coupling procedure with planar interpolation operators between adjoining nonconforming elements is used. Curvilinear volume metric terms are numerically approximated via a minimization procedure and satisfy the discrete geometric conservation law conditions. Distinct curvilinear surface metrics are used on the adjoining interfaces to construct the interface coupling terms, thereby localizing the discrete geometric conservation law constraints to each individual element. The resulting scheme is entropy conservative/stable, element-wise conservative, and freestream preserving. Viscous interface dissipation operators are developed that retain the entropy stability of the base scheme. The accuracy and stability properties of the resulting numerical scheme are shown to be comparable to those of the original conforming scheme (achieving ~p+1 convergence) in the context of the viscous shock problem, the Taylor-Green vortex problem at a Reynolds number of Re=1,600, and a subsonic turbulent flow past a sphere at Re = 2,000.
We introduce solution dependent finite difference stencils whose coefficients adapt to the current numerical solution by minimizing the truncation error in the least squares sense. The resulting scheme has the resolution capacity of dispersion relation preserving difference stencils in under-resolved domains, together with the high order convergence rate of conventional central difference methods in well resolved regions. Numerical experiments reveal that the new stencils outperform their conventional counterparts on all grid resolutions from very coarse to very fine.
This work examines the development of an entropy conservative (for smooth solutions) or entropy stable (for discontinuous solutions) space–time discontinuous Galerkin (DG) method for systems of nonlinear hyperbolic conservation laws. The resulting numerical scheme is fully discrete and provides a bound on the mathematical entropy at any time according to its initial condition and boundary conditions. The crux of the method is that discrete derivative approximations in space and time are summation-by-parts (SBP) operators. This allows the discrete method to mimic results from the continuous entropy analysis and ensures that the complete numerical scheme obeys the second law of thermodynamics. Importantly, the novel method described herein does not assume any exactness of quadrature in the variational forms that naturally arise in the context of DG methods. Typically, the development of entropy stable schemes is done on the semidiscrete level ignoring the temporal dependence. In this work, we demonstrate that creating an entropy stable DG method in time is similar to the spatial discrete entropy analysis, but there are important (and subtle) differences. Therefore, we highlight the temporal entropy analysis throughout this work. For the compressible Euler equations, the preservation of kinetic energy is of interest besides entropy stability. The construction of kinetic energy preserving (KEP) schemes is, again, typically done on the semidiscrete level similar to the construction of entropy stable schemes. We present a generalization of the KEP condition from Jameson to the space–time framework and provide the temporal components for both entropy stability and kinetic energy preservation. The properties of the space–time DG method derived herein are validated through numerical tests for the compressible Euler equations. Additionally, we provide, in appendices, how to construct the temporal entropy stable components for the shallow water or ideal magnetohydrodynamic (MHD) equations.
The entropy conservative/stable staggered grid tensor-product algorithm of Parsani et al. [1] is extended to multidimensional SBP discretizations. The required SBP preserving interpolation operators are proven to exist under mild restrictions and the resulting algorithm is proven to be entropy conservative/stable as well as elementwise conservative. For 2-dimensional simplex elements, the staggered grid algorithm is shown to be more accurate and have a larger maximum time step restriction as compared to the collocated algorithm. The staggered algorithm significantly reduces the number of (computationally expensive) two-point flux evaluations, which is potentially important for both explicit and implicit time-marching schemes. Furthermore, the staggered algorithm requires fewer degrees of freedom for comparable accuracy, which has favorable implications for implicit time-marching schemes.
We consider the numerical simulation of the acoustic wave equations arising from seismic applications, for which staggered grid finite difference methods are popular choices due to their simplicity and efficiency. We relax the uniform grid restriction on finite difference methods and allow the grids to be block-wise uniform with nonconforming interfaces. In doing so, variations in the wave speeds of the subterranean media can be accounted for more efficiently. Staggered grid finite difference operators satisfying the summation-by-parts (SBP) property are devised to approximate the spatial derivatives appearing in the acoustic wave equation. These operators are applied within each block independently. The coupling between blocks is achieved through simultaneous approximation terms (SATs), which impose the interface conditions weakly, i.e., by penalty. Ratio of the grid spacing of neighboring blocks is allowed to be rational number, for which specially designed interpolation formulas are presented. These interpolation formulas constitute key pieces of the simultaneous approximation terms. The overall discretization is shown to be energy-conserving and examined on test cases of both theoretical and practical interests, delivering accurate and stable simulation results.
New entropy stable spectral collocation schemes of arbitrary order of accuracy are developed for the unsteady 3-D Euler and Navier-Stokes equations on dynamic unstructured grids. To take into account the grid motion and deformation, we use an arbitrary Lagrangian-Eulerian formulation. As a result, moving and deforming hexahedral grid elements are individually mapped onto a cube in the fixed reference system of coordinates. The proposed scheme is constructed by using the skew-symmetric form of the Navier-Stokes equations, which are discretized by using summation-by-parts spectral collocation operators that preserve the conservation properties of the original governing equations. Furthermore, the metric coefficients are approximated such that the geometric conservation laws are satisfied exactly on both static and dynamic grids. To make the scheme entropy stable, a new entropy conservative flux is derived for the 3-D Euler and Navier-Stokes equations on dynamic unstructured grids. The new flux preserves the design order of accuracy of the original spectral collocation scheme and guarantees entropy conservation on moving and deforming grids. We present numerical results demonstrating design order of accuracy and freestream preservation properties of the new schemes for both the Euler and Navier-Stokes equations on dynamic unstructured grids.
The construction of high order entropy stable collocation schemes on quadrilateral and hexahedral elements has relied on the use of Gauss-Legendre-Lobatto collocation points and their equivalence with summation-by-parts (SBP) finite difference operators. In this work, we show how to efficiently generalize the construction of semi-discretely entropy stable schemes on tensor product elements to Gauss points and generalized SBP operators. Numerical experiments suggest that the use of Gauss points significantly improves accuracy on curved meshes.
Methodologies are presented that enable the construction of provably linearly stable and conservative high-order discretizations of partial differential equations in curvilinear coordinates based on generalized summation-by-parts operators, including operators with dense-norm matrices. Specifically, three approaches are presented for the construction of stable and conservative schemes in curvilinear coordinates using summation-by-parts (SBP) operators that have a diagonal norm but may or may not include boundary nodes: (1) the mortar-element approach, (2) the global SBP-operator approach, and (3) the staggered-grid approach. Moreover, the staggered-grid approach is extended to enable the development of stable dense-norm operators in curvilinear coordinates. In addition, collocated upwind simultaneous approximation terms for the weak imposition of boundary conditions or inter-element coupling are extended to curvilinear coordinates with the new approaches. While the emphasis in the paper is on tensor-product SBP operators, the approaches that are covered are directly applicable to multidimensional SBP operators.
We present an entropy stable numerical scheme subject to no-slip wall boundary conditions. To enforce entropy stability only the no-penetration boundary condition and a temperature condition are needed at a wall, and this leads to an L-2 bound on the conservative variables. In this article, we take the next step and design a finite difference scheme that also bounds the velocity gradients. This necessitates the use of the full no-slip conditions.
We present and analyze an entropy-stable semi-discretization of the Euler equations based on high-order summation-by-parts (SBP) operators. In particular, we consider general multidimensional SBP elements, building on and generalizing previous work with tensor–product discretizations. In the absence of dissipation, we prove that the semi-discrete scheme conserves entropy; significantly, this proof of nonlinear L2 stability does not rely on integral exactness. Furthermore, interior penalties can be incorporated into the discretization to ensure that the total (mathematical) entropy decreases monotonically, producing an entropy-stable scheme. SBP discretizations with curved elements remain accurate, conservative, and entropy stable provided the mapping Jacobian satisfies the discrete metric invariants; polynomial mappings at most one degree higher than the SBP operators automatically satisfy the metric invariants in two dimensions. In three-dimensions, we describe an elementwise optimization that leads to suitable Jacobians in the case of polynomial mappings. The properties of the semi-discrete scheme are verified and investigated using numerical experiments.
Two new implicit–explicit, additive Runge–Kutta (ARK2) methods are given with fourth- and fifth-order formal accuracies, respectively. Both combine explicit Runge–Kutta (ERK) methods with explicit, singly-diagonally implicit Runge–Kutta (ESDIRK) methods and include an embedded method for error control. The two methods have ESDIRKs which are internally L-stable on stages three and higher, have only modestly negative eigenvalues to the stage and step algebraic-stability matrices and have stage-order two. To improve computational efficiency, the fourth-order method has a diagonal coefficient of 0.1235. This is done to offset much of the extra computational cost of an extra stage by facilitating iterative convergence at each stage. Linear stability domains for both ERK methods have been made quite large and the dominant coupling stability term between the stability of the ESDIRK and ERK for very stiff modes has been removed. Though the fourth-order method is one of the best all-around fourth-order IMEX ARK2 of which we are aware, the fifth-order method is likely best suited to mildly stiff problems with tight error tolerances. Methods are tested using the Van der Pol and Kaps' singular-perturbation problems. Results suggest that these new methods represent an improvement over existing methods of the same class.
This work presents an entropy stable discontinuous Galerkin (DG) spectral element approximation for systems of non-linear conservation laws with general geometric (h) and polynomial order (p) non-conforming rectangular meshes. The crux of the proofs presented is that the nodal DG method is constructed with the collocated Legendre-Gauss-Lobatto nodes. This choice ensures that the derivative/mass matrix pair is a summation-by-parts (SBP) operator such that entropy stability proofs from the continuous analysis are discretely mimicked. Special attention is given to the coupling between nonconforming elements as we demonstrate that the standard mortar approach for DG methods does not guarantee entropy stability for non-linear problems, which can lead to instabilities. As such, we describe a precise procedure and modify the mortar method to guarantee entropy stability for general non-linear hyperbolic systems on h/p non-conforming meshes. We verify the high-order accuracy and the entropy conservation/stability of fully non-conforming approximation with numerical examples.
We present an entropy-stable semi-discretization of the Euler equations. The scheme is based on high-order summation-by-parts (SBP) operators for triangular and tetrahedral elements, although the theory is applicable to multidimensional SBP operators on more general elements. While there are established methods for proving stability of linear equations, such as energy analysis, they are not adequate for nonlinear equations. To address nonlinear stability, we use the matrix properties of the SBP operators combined with entropy-conserving numerical flux functions. This allows us to prove that the semi-discrete scheme conserves entropy. Significantly, the proof does not rely on integral exactness, and, therefore, the discretization has a stronger claim of robustness than a similar finite-element method. The addition of an upwinded term to the entropy-conservative scheme makes it entropy-stable. This generalizes previous work proving entropy stability for tensor-product elements to more general elements, including simplex elements. Numerical experiments are conducted to verify accuracy and entropy conservation on an isentropic vortex flow.
High-order numerical methods that satisfy a discrete analog of the entropy inequality are uncommon. Indeed, no proofs of nonlinear entropy stability currently exist for high-order weighted essentially nonoscillatory (WENO) finite volume or weak-form finite element methods. Herein, a new family of fourth-order WENO spectral collocation schemes is developed, that are nonlinearly entropy stable for the one-dimensional compressible Navier–Stokes equations. Individual spectral elements are coupled using penalty type interface conditions. The resulting entropy stable WENO spectral collocation scheme achieves design order accuracy, maintains the WENO stencil biasing properties across element interfaces, and satisfies the summation-by-parts (SBP) operator convention, thereby ensuring nonlinear entropy stability in a diagonal norm. Numerical results demonstrating accuracy and nonoscillatory properties of the new scheme are presented for the one-dimensional Euler and Navier–Stokes equations for both continuous and discontinuous compressible flows.