High-order Shock-capturing Schemes Research Articles

A novel finite-difference converged ENO scheme for steady-state simulations of Euler equations

The high-order shock-capturing scheme is critical for compressible fluid simulations, in particular for cases where strong shockwaves are present. However, for the steady-state solution of Euler equations, the high-order shock-capturing schemes usually suffer from the difficulty that the numerical residual hangs at a relatively high level instead of converging to machine zero. In this paper, a new paradigm of finite-difference converged ENO (CENO) scheme for steady-state simulations of Euler equations is proposed. Different from the existing methods dedicating to eliminating the slight post-shock oscillation, the unsteady effect of shock-induced numerical instability on the nonlinear weighting process as well as the resultant variation of the numerical formula for the spatially fixed cell interface flux at different time instants are demonstrated to be the core contributors to the non-convergence of high-order shock-capturing schemes. Unlike classical ENO/WENO/TENO schemes, which generate the final reconstruction based only on spatial information, evolutionary information is also exploited by a novel converged stencil selection strategy in the present scheme. With this strategy, the unsteady effect of shock-induced numerical instability on the nonlinear weighting process is incorporated and eliminated by a tailored adaptive scale-separation algorithm, while shock-capturing capabilities are not sacrificed. A set of well-established and newly designed challenging benchmark simulations involving strong shockwaves, contact discontinuities and rarefaction waves demonstrates that the new fifth-order CENO5 scheme features superior convergence properties and achieves the essentially non-oscillation property better when compared to the robust TENO5 scheme and the state-of-the-art WENO-type schemes designed for steady-state problems.

Steady-state simulation of Euler equations by the discontinuous Galerkin method with the hybrid limiter

In the realm of steady-state solutions of Euler equations, the challenge of achieving convergence of residue close to machine zero has long plagued high-order shock-capturing schemes, particularly in the presence of strong shock waves. In response to this issue, this paper presents a hybrid limiter designed specifically for the discontinuous Galerkin (DG) method which incorporates a pseudo-time marching method. This hybrid limiter, initially introduced in DG methods for unsteady problems, preserves the local data structure while enhancing resolution through effective and precise shock capturing. Notably, for the steady problem, the hybridization of the DG solution with the cell average is employed, deviating from the low-order limited DG solution typically used for unsteady problems. This approach seamlessly integrates the previously distinct troubled cell indicator and limiter components resulting in a more cohesive and efficient limiter. It obviates the need for characteristic decomposition and intercell communication, leading to substantial reductions in computational costs and enhancement in parallel efficiency. Numerical experiments are presented to demonstrate the robust performance of the hybrid limiter for steady-state Euler equations both on the structured and unstructured meshes.

A new high-order shock-capturing TENO scheme combined with skew-symmetric-splitting method for compressible gas dynamics and turbulence simulation

The high-order shock-capturing scheme is one of the main building blocks for the simulation of the compressible fluid characterized by strong shockwaves and broadband length scales. However, the classical shock-capturing scheme fails to perform long-time stable and non-dissipative simulations since the quadratic invariants associated with the conservation equations cannot be conserved as a result of the inherent numerical dissipation. Additionally, the overall computational cost for classical shock-capturing schemes is quite expensive as a result of the time-consuming local characteristic decomposition and the nonlinear-weights computing process. In this work, based on a new efficient discontinuity indicator, which distinguishes the non-smooth high-wavenumber fluctuations and discontinuities from smooth scales in the wavenumber space, a paradigm of high-order shock-capturing scheme by recasting the non-dissipative skew-symmetric-splitting method with newly optimized dispersion property for smooth flow scales and invoking the nonlinear targeted ENO (TENO) schemes for non-smooth discontinuities is proposed. The resulting TENO-S scheme not only successfully performs long-time stable computations for smooth flows without numerical dissipation, but also recovers the robust shock-capturing capabilities with adaptive numerical dissipation. Without the necessity of parameter tuning case by case, extensive benchmark simulations involving a wide range of flow length scales and strong discontinuities demonstrate that the proposed TENO-S scheme performs significantly better than the straightforward deployment of WENO/TENO-family schemes with better spectral property and higher computational efficiency.

BROADCAST: A high-order compressible CFD toolbox for stability and sensitivity using Algorithmic Differentiation

The evolution of any complex dynamical system is described by its state derivative operators. However, the extraction of the exact N-order state derivative operators is often inaccurate and requires approximations. The open-source CFD (Computational Fluid Dynamics) code called BROADCAST discretises the compressible Navier-Stokes equations and then extracts the linearised N-derivative operators through Algorithmic Differentiation (AD) providing a toolbox for laminar flow dynamics. Furthermore, the gradients through adjoint derivation are extracted either by transposition of the linearised operator or through the backward mode of the AD tool. The software includes base-flow computation and linear global stability analysis via eigen-decomposition of the linearised operator or via singular value decomposition of the resolvent operator. Sensitivity tools as well as weakly nonlinear analysis complete the package. The numerical method for the spatial discretisation of the equations consists of a finite-difference high-order shock-capturing scheme applied within a finite volume framework on 2D curvilinear structured grids. The stability and sensitivity tools are demonstrated on two cases: a cylinder flow at low Mach number and a hypersonic boundary layer. Program summaryProgram Title: BROADCASTCPC Library link to program files:https://doi.org/10.17632/hcf893cr2n.1Developer's repository link:https://github.com/OneraHub/BroadcastLicensing provisions: Mozilla Public License 2.0Programming language: Python, FortranNature of problem: The extraction of the exact N-order state derivative operators useful to study complex dynamical systems is often inaccurate and requires approximations from the user.Solution method: This program discretises the compressible Navier-Stokes equations and then extracts the linearised state derivative operators through Algorithmic Differentiation providing a toolbox for laminar flow dynamics.

Symmetry-preserving enforcement of low-dissipation method based on boundary variation diminishing principle

A class of high-order shock-capturing schemes, PnTm-BVD (Deng et al., 2019; Deng et al., 2020) schemes, have been devised to solve the Euler equations with substantially reduced numerical dissipation, which enable high-resolution simulations to resolve flow structures of wider range scales. In such simulations with low dissipation, errors of round-off level might grow and contaminate the numerical solutions. A typical example of such problems is the loss of symmetry in the numerical solutions for physical problems of symmetric configurations even if the schemes are mathematically in line with the symmetry rules. In this study, the mechanisms of symmetry-breaking in a finite volume framework with the P4T2-BVD reconstruction scheme are thoroughly examined. Particular attention has been paid to remove the possible causes due to the lack of associativity in floating-point arithmetic which is associated with round-off errors. Modifications and new techniques are proposed to completely remove the possible causes for symmetry breaking in different components of the P4T2-BVD finite volume solver. Benchmark tests that have symmetric solution structures are used to verify the proposed methods. The numerical results demonstrate the perfect symmetric solution structures.

High-order numerical scheme for compressible multi-component real gas flows using an extension of the Roe approximate Riemann solver and specific Monotonicity-Preserving constraints

The purpose of this paper is to develop a high-order shock-capturing scheme capable of predicting flows where shock waves with high-temperature jumps interact with multi-component real gas mixtures, assuming a local thermodynamic equilibrium. We first propose a generalization of the Roe solver for distinct species with non-ideal thermodynamic properties that relies on the original method proposed by Vinokur & Montagné [1]. This method uses an approximation of compressibility factors to estimate a coherent value of the speed of sound at the Roe averaged state.This Roe averaged state is introduced in the One-Step Monotonicity-Preserving (OSMP) scheme, originally developed by Daru and Tenaud [2], to obtain an extension to the high-order with Lax-Wendroff procedure adequate for dealing with non-ideal gas flows. To avoid thermodynamic inconsistencies in the evolution of the Roe average state over a large stencil, we propose to reformulate the discrete total energy flux of the initial solver. This new formulation uses a combination of Riemann invariants related to the species mass fractions and avoids the influence of the independent values of the compressibility factors in the total energy flux computation. An additional M-P constraint on this new combination allows dealing with discontinuities. Based on the averaged speed of sound estimated by our proposed extension of the Vinokur & Montagné method, we demonstrate that this new formulation is equivalent to selecting a new combination of compressibility factors that completely fulfill the jump relationships of the Riemann problem.To properly capture discontinuities while optimizing the number of numerical cells, the new high-order OSMP scheme is combined with an Adaptive Multiresolution [3] procedure to automatically refine grid in regions where steep gradients occur and coarsen grid elsewhere. The order of the numerical method is evaluated on the convection of density and mass fraction waves. Its capability of capturing discontinuities is validated on a 1-D shock tube problem with a mixture of Nitrogen, Oxygen and dense refrigerant R22 gases. We show that smooth solutions, as well as discontinuities, are recovered with high accuracy. The 2-D interaction between a shock wave in Air with a cylindrical bubble initially filled with dense refrigerant R22 gas is also considered. Present results compare very well with both a recent fully resolved numerical solution of ideal gases and experimental results obtained with real gases. Compared to ideal gas solutions corresponding to calorically perfect gas, drastic changes are recorded on the predicted temperature and the bubble flow patterns that fully justify the use of relevant thermodynamics and the proposed numerical method to account for real gas properties.

Numerical symmetry-preserving techniques for low-dissipation shock-capturing schemes

Modern applications of computational fluid dynamics involve complex interactions across scales such as shock interactions with turbulent structures and multiphase interfaces. Such phenomena, which occur at very small physical viscosity, require high-resolution and low-dissipation compressible flow solvers. Many recent publications have focused on the design of high-order accurate numerical schemes and provide e.g. weighted essentially non-oscillatory (WENO) stencils up to 17th order for this purpose. As shown in detail by different authors, such schemes tremendously decrease adverse effects of numerical dissipation. However, such schemes are prone to numerically induced symmetry breaking which renders validation for the targeted problem range problematic.In this paper, we show that symmetry-breaking relates to vanishing numerical viscosity and is driven systematically by algorithmic floating-point effects which are no longer hidden by numerical dissipation. We propose a systematic procedure to deal with such errors by numerical and algorithmic formulations which respect floating-point arithmetic. We show that by these procedures inherent symmetries are preserved for a broad range of test cases with high-order shock-capturing schemes in particular in the high-resolution limit for both 2D and 3D.

Analysis of mild ignition in a shock tube using a highly resolved 3D-LES and high-order shock-capturing schemes

A highly resolved three-dimensional large-eddy simulation (LES) is presented for a shock tube containing a stoichiometric hydrogen–oxygen ( $$\hbox {H}_2$$ / $$\hbox {O}_2$$ ) mixture, and the results are compared against experimental results. A parametric study is conducted to test the effects of grid resolution, numerical scheme, and initial conditions before the 3D simulations are presented in detail. An approximate Riemann solver and a high-order interpolation scheme are used to solve the conservation equations of the viscous, compressible fluid and to account for turbulence behind the reflected shock. Chemical source terms are calculated by a finite-rate model. Simultaneous results of pseudo-Schlieren, temperature, pressure, and species are presented. The ignition delay time is predicted in agreement with the experiments by the three-dimensional simulations. The mechanism of mild ignition is analysed by Lagrangian tracer particles, tracking temperature histories of material particles. We observed strongly increased temperatures in the core region away from the end wall, explaining the very early occurrence of mild ignition in this case.

Modeling Bed Evolution Using Weakly Coupled Phase-Resolving Wave Model and Wave-Averaged Sediment Transport Model

In this paper, we propose a model for the simulation of the bed evolution dynamics in coastal regions characterized by articulated morphologies. An integral form of the fully nonlinear Boussinesq equations in contravariant formulation, in which Christoffel symbols are absent, is proposed in order to simulate hydrodynamic fields from deep water up to just seaward of the surf zones. Breaking wave propagation in the surf zone is simulated by integrating the nonlinear shallow water equations with a high-order shock-capturing scheme. The near-bed instantaneous flow velocity and the intra-wave hydrodynamic quantities are calculated by the momentum equation integrated over the turbulent boundary layer. The bed evolution dynamics is calculated starting from the contravariant formulation of the advection-diffusion equation for the suspended sediment concentration in which the advective sediment transport terms are formulated according to a quasi-three-dimensional approach, and taking into account the contribution given by the spatial variation of the bed load transport. The model is validated against several tests by comparing numerical results with experimental data. The ability of the proposed model to represent the sediment transport phenomena in a morphologically articulated coastal region is verified by numerically simulating the long-term bed evolution in the coastal region opposite Pescara harbor (in Italy) and comparing numerical results with the field data.

Numerical dissipation control in high order shock-capturing schemes for LES of low speed flows

The Yee & Sjögreen adaptive numerical dissipation control in high order scheme (High Order Filter Methods for Wide Range of Compressible Flow Speeds, ICOSAHOM 09, 2009) is further improved for DNS and LES of shock-free turbulence and low speed turbulence with shocklets. There are vastly different requirements in the minimization of numerical dissipation for accurate turbulence simulations of different compressible flow types and flow speeds. Traditionally, the method of choice for shock-free turbulence and low speed turbulence are by spectral, high order central or high order compact schemes with high order linear filters. With a proper control of a local flow sensor, appropriate amount of numerical dissipation in high order shock-capturing schemes can have spectral-like accuracy for compressible low speed turbulent flows. The development of the method includes an adaptive flow sensor with automatic selection on the amount of numerical dissipation needed at each flow location for more accurate DNS and LES simulations with less tuning of parameters for flows with a wide range of flow speed regime during the time-accurate evolution, e.g., time varying random forcing. An automatic selection of the different flow sensors catered to the different flow types is constructed. A Mach curve and high-frequency oscillation indicators are used to reduce the tuning of parameters in controlling the amount of shock-capturing numerical dissipation to be employed for shock-free turbulence, low speed turbulence and turbulence with strong shocks. In Kotov et al. (High Order Numerical Methods for LES of Turbulent Flows with Shocks, ICCFD8, Chengdu, Sichuan, China, July 14–18, 2014) the LES of a turbulent flow with a strong shock by the Yee & Sjögreen scheme indicated a good agreement with the filtered DNS data. A work in progress for the application of the adaptive flow sensor for compressible turbulence with time-varying random forcing is forthcoming. The present study examines the versatility of the Yee & Sjögreen scheme for DNS and LES of traditional low speed flows without forcing. Special attention is focused on the accuracy performance of this scheme using the Smagorinsky and the Germano–Lilly SGS models.

Hybrid spectral difference/embedded finite volume method for conservation laws

Recently, interests have been increasing towards applying the high-order methods to various engineering applications with complex geometries [30]. As a result, a family of discontinuous high-order methods, such as Discontinuous Galerkin (DG), Spectral Volume (SV) and Spectral Difference (SD) methods, is under active development. These methods provide high-order accurate solutions and are highly parallelizable due to the local solution reconstruction within each element. But, these methods suffer from the Gibbs phenomena when discontinuities are present in the flow fields. Various types of limiters [43–45] and artificial viscosity [46,48] have been employed to overcome this problem.A novel hybrid spectral difference/embedded finite volume method is introduced in order to apply a discontinuous high-order method for large scale engineering applications involving discontinuities in the flows with complex geometries. In the proposed hybrid approach, the finite volume (FV) element, consisting of structured FV subcells, is embedded in the base hexahedral element containing discontinuity, and an FV based high-order shock-capturing scheme is employed to overcome the Gibbs phenomena. Thus, a discontinuity is captured at the resolution of FV subcells within an embedded FV element. In the smooth flow region, the SD element is used in the base hexahedral element. Then, the governing equations are solved by the SD method. The SD method is chosen for its low numerical dissipation and computational efficiency preserving high-order accurate solutions. The coupling between the SD element and the FV element is achieved by the globally conserved mortar method [56]. In this paper, the 5th-order WENO scheme with the characteristic decomposition is employed as the shock-capturing scheme in the embedded FV element, and the 5th-order SD method is used in the smooth flow field.The order of accuracy study and various 1D and 2D test cases are carried out, which involve the discontinuities and vortex flows. Overall, it is shown that the proposed hybrid method results in comparable or better simulation results compared with the standalone WENO scheme when the same number of solution DOF is considered in both SD and FV elements.

Implementation of a Robust Weighted Compact Nonlinear Scheme for Modeling of Hydrogen/Air Detonation

In order to simulate detonation, a high-order shock-capturing scheme that models chemical reactions is implemented, and its resolution is examined by testing it on several numerical problems. A robust weighted compact nonlinear scheme (RWCNS) is adopted to take advantage of its robustness and ability to handle different flux types. This study shows the high resolution of the RWCNS compared with the conventional scheme. The results show that the RWCNS predicts the detailed vortex structure behind the detonation wavefront.

A Comparative Study of Accuracy of Shock Capturing Schemes for Simulation of Shock/Acoustic Wave Interactions

We simulate transmission of a small-amplitude disturbance wave through a shock wave. Results of our numerical experiments performed with different high-order shock-capturing schemes show that the capability of a scheme to correctly predict the amplification of the disturbances crucially depends on the Riemann solver used in evaluation of numerical fluxes. Incorrectly high amplification rates are produced by the solvers resolving shock waves sharply, with no interior points in numerical profiles of steady shock waves. In particular, both the exact Riemann solver and the popular Roe flux difference splitting demonstrate such unphysical behavior. A possible explanation of such behavior is proposed. More dissipative solvers, such as the global Lax–Friedrichs splitting, produce transmission coefficients close to the predictions of linear theory.

Computational challenges for simulations related to the NASA electric arc shock tube (EAST) experiments

The goal of this study is to gain some physical insights and an understanding of the computational challenges for the simulations related to the hypersonic nonequilibrium multi-species and multi-reaction experiments on the NASA Electric Arc Shock Tube (EAST). While experimental measurement does not provide any information about the radial structure of this type of flow, accurate and reliable numerical simulations can provide more insight into the physical structure of the flow to aid the design of atmospheric entry spacecrafts. The paper focuses on the spurious numerics which take place in numerical simulations of the subject physics containing stiff source terms and discontinuities. This paper is based on the knowledge gained from Yee et al. on simple reacting test cases (Yee et al. 2013, [9]) as a guide to reveal the computational challenges involved for such an extreme flow type. The results of the 1D and 2D EAST viscous and inviscid simulations using a simplified physical model are presented. The computation reveals, for the first time, that the 2D viscous model which contains both shocks and shears exhibits Tollmien–Schlichting-like instability complex patterns at the boundary layer. In addition to exhibiting spurious numerical behavior of wrong propagation speed of discontinuities by typical shock-capturing methods, there is improved understanding on the cause of numerical difficulties by previous investigators. One example is that the relative distance between the shocks and shear/contact is different from one grid spacing to another for each considered high order shock-capturing scheme. The results presented can provide insight on the numerical instability observed by previous investigations and future algorithm development for this type of extreme flow.

Spurious behavior of shock-capturing methods by the fractional step approach: Problems containing stiff source terms and discontinuities

The goal of this paper is to relate numerical dissipations that are inherited in high order shock-capturing schemes with the onset of wrong propagation speed of discontinuities. For pointwise evaluation of the source term, previous studies indicated that the phenomenon of wrong propagation speed of discontinuities is connected with the smearing of the discontinuity caused by the discretization of the advection term. The present study focuses only on solving the reactive system by the fractional step method using the Strang splitting. Studies shows that the degree of wrong propagation speed of discontinuities is highly dependent on the accuracy of the numerical method. The manner in which the smearing of discontinuities is contained by the numerical method and the overall amount of numerical dissipation being employed play major roles. Depending on the numerical method, time step and grid spacing, the numerical simulation may lead to (a) the correct solution (within the truncation error of the scheme), (b) a divergent solution, (c) a wrong propagation speed of discontinuities solution or (d) other spurious solutions that are solutions of the discretized counterparts but are not solutions of the governing equations. The findings might shed some light on the reported difficulties in numerical combustion and problems with stiff nonlinear (homogeneous) source terms and discontinuities in general.

Comparative Study of Three High Order Schemes for LES of Temporally Evolving Mixing Layers

AbstractThree high order shock-capturing schemes are compared for large eddy simulations (LES) of temporally evolving mixing layers for different convective Mach numbers ranging from the quasi-incompressible regime to highly compressible supersonic regime. The considered high order schemes are fifth-order WENO (WENO5), seventh-order WENO (WENO7) and the associated eighth-order central spatial base scheme with the dissipative portion of WENO7 as a nonlinear post-processing filter step (WENO7fi). This high order nonlinear filter method of Yee & Sjögreen is designed for accurate and efficient simulations of shock-free compressible turbulence, turbulence with shocklets and turbulence with strong shocks with minimum tuning of scheme parameters. The LES results by WENO7fi using the same scheme parameter agree well with experimental results compiled by Barone et al., and published direct numerical simulations (DNS) work of Rogers & Moser and Pantano & Sarkar, whereas results by WENO5 and WENO7 compare poorly with experimental data and DNS computations.

Long-time simulations of the Kelvin-Helmholtz instability using an adaptive vortex method

The nonlinear evolution of an interface subject to a parallel shear flow is studied by the vortex sheet model. We perform long-time computations for the vortex sheet in density-stratified fluids by using the point vortex method and investigate late-time dynamics of the Kelvin-Helmholtz instability. We apply an adaptive point insertion procedure and a high-order shock-capturing scheme to the vortex method to handle the nonuniform distribution of point vortices and enhance the resolution. Our adaptive vortex method successfully simulates chaotically distorted interfaces of the Kelvin-Helmholtz instability with fine resolutions. The numerical results show that the Kelvin-Helmholtz instability evolves a secondary instability at a late time, distorting the internal rollup, and eventually develops to a disordered structure.

Wave processes in the shock layer on a flat plate at an angle of attack

A numerical and experimental study of receptivity of the viscous shock layer on a flat plate aligned at an angle of attack to external acoustic perturbations is performed. Density and pressure fluctuations are measured in experiments at the free-stream Mach number M ∞ = 21 and Reynolds number Re 1 = 6·10 5 m −1 . Direct numerical simulations of receptivity of the viscous shock layer to external acoustic perturbations in wide ranges of the governing parameters are performed by solving the Navier-Stokes equations with the use of high-order shock-capturing schemes. The calculated intensities of density and pressure fluctuations are found to be in good agreement with experimental data. Results of the study show that entropy-vortex disturbances dominate in the shock layer at small angles of attack, whereas acoustic perturbations prevail at angles of attack above 20°.

Wave processes in a viscous shock layer and control of fluctuations

Generation and development of disturbances in a hypersonic viscous shock layer on a flat plate is studied both experimentally and numerically. The study is performed at the Mach number M∞ = 21 and the Reynolds number ReL = 1.44 × 105 and is aimed at elucidating the physical mechanisms that govern the receptivity and instability of the shock layer at extremely high hypersonic velocities. The experiments are conducted in a hypersonic nitrogen-driven wind tunnel. An electron-beam fluorescence technique, a Pitot probe and a piezoceramic transducer are used to measure the mean density and Mach number contours, as well as density and pressure fluctuations, their spectra and spatial distributions in the shock layer. Direct numerical simulations are performed by solving the Navier–Stokes equations with a high-order shock-capturing scheme in a computational domain including the leading and trailing edges of the plate, so that the bow shock wave and the wake behind the plate are also simulated. It is demonstrated that computational and experimental data characterizing the mean flow field, intensity of density fluctuations and their spatial distributions in the shock layer are in close agreement. It is found that excitation of the shock layer by external acoustic waves leads to generation of entropy–vortex disturbances with two maxima of density fluctuations: directly behind the shock wave and on the external edge of the boundary layer. At the same time, the pressure fluctuations decay inward into the shock layer, away from the shock, which agrees with the linear theory of interaction of shock waves with small perturbations. Thus, the entropy–vortex disturbances are shown to dominate in the hypersonic shock layer at very high Mach numbers, in contrast with the boundary layers at moderate hypersonic velocities where acoustic modes are most important. A parametric numerical study of wave processes in the shock layer induced by external acoustic waves is performed with variations of frequency, amplitude and angle of propagation of external disturbances. The amplitude of generated disturbances is observed to grow and decay periodically along the streamwise coordinate, and the characteristics of these variations depend on the frequency and direction of incident acoustic waves. The hypersonic shock layer excited by periodic blowing and suction near the leading edge is also investigated; in the experiments, this type of excitation is obtained by using an oblique-cut whistle. It is shown that blowing/suction generates disturbances resembling those generated by external acoustic waves, with similar spatial distributions and phase velocities. This result paves the way for active control of instability development in the shock layer by means of destructive interference of two types of disturbances. Numerical simulations are performed to show that instability waves can be significantly amplified or almost entirely suppressed, depending on the relative phase of blowing/suction and acoustic disturbances. Wind-tunnel experiments completely confirm this numerical prediction. Thus, the feasibility of delaying instability development in the hypersonic shock layer has been demonstrated for the first time.

Localized artificial diffusivity scheme for discontinuity capturing on curvilinear meshes

A simple and efficient localized artificial diffusivity scheme is developed for the purpose of capturing discontinuities on curvilinear and anisotropic meshes using a high-order compact differencing scheme. The artificial diffusivity is dynamically localized in space to capture different types of discontinuities such as a shock wave, contact surface or material discontinuity. The method is intended for use with large-eddy simulation of compressible transitional and turbulent flows. The method captures the discontinuities on curvilinear and anisotropic meshes with minimum impact on the smooth flow regions. The amplitude of wiggles near a discontinuity and the number of grid points used to capture the discontinuity do not depend on the mesh size. The comparisons between the proposed method and high-order shock-capturing schemes illustrate the advantage of the method for the simulation of flows involving shocks, turbulence and their interactions. The multi-dimensional formulation is tested on a variety of 1D and 2D, steady and unsteady, different types of discontinuity-related problems on curvilinear and anisotropic meshes. A simplification of the method which reduces the computational cost does not show any major detrimental effect on the discontinuity capturing under the conditions examined.